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I.  Introduction 

We  are  in  the  second  year  of  a  three  year  effort  focused  on  developing  enhanced 
contrast  thermal  acoustic  imaging  (TAI)  technology  for  the  detection  of  breast  cancer  by 
combining  amplitude-modulated  (AM)  electromagnetic  (EM)  field  excitation,  resonant 
acoustic  scattering,  and  advanced  signal  processing  techniques.  EM-induced  TAI 
combines  the  merits  of  both  EM  stimulation  and  ultrasound  imaging,  while  overcoming 
their  respective  limitations.  EM  imaging  provides  excellent  contrast  between  cancerous 
and  normal  breast  tissue,  but  the  long  wavelengths  provide  poor  spatial  resolution. 
Conventional  ultrasound  imaging  possesses  very  fine  millimeter-range  spatial  resolution 
but  poor  soft  tissue  contrast.  While  EM-induced  TAI  possesses  great  promise,  the  thermal 
acoustic  signals  tend  to  be  weak.  However,  when  the  tumor  is  excited  into  resonance  via 
EM  stimulation,  the  effective  acoustic  scattering  cross-section  may  increase  by  a  factor  in 
excess  of  100  based  on  predictions  for  microsphere-based  ultrasound  contrast  agents. 
Such  an  increase  would  truly  be  revolutionary,  making  the  EM-induced  TAI  technology  a 
very  promising  candidate  for  routine  breast  screening.  To  induce  the  resonant  response 
from  the  tumor,  we  are  considering  various  approaches  including,  for  example,  AM 
continuous  wave  (CW)  EM  stimulation,  where  the  modulation  frequency  range  contains 
the  predicted  resonant  frequencies  for  a  distribution  of  tumor  sizes  and  contrast  ratios. 
The  carrier  frequency  of  the  EM  stimulation  can  be  fixed  and  chosen  for  the  best 
penetration  and  heat  absorption.  The  image  formation  methods  in  the  existing  TAI 
systems  are  predominantly  data-independent  delay-and-sum  (with  or  without  weighting) 
type  of  approaches.  These  approaches  tend  to  have  poor  resolution  (relative  to  the  best 
possible  resolution  a  transducer  array  can  offer)  and  high  sidelobe  problems,  especially 
when  the  transducer  array  is  not  composed  of  uniformly  and  linearly  spaced  transducers, 
which  is  the  case  for  the  existing  TAI  systems.  We  are  devising  adaptive  image  formation 
algorithms  to  achieve  high  resolution  and  excellent  interference  and  noise  suppression 
capability. 

The  research  efforts  this  year  were  focused  on  the  following  areas:  initial  phantom 
resonance  experiments,  further  development  of  tissue  mimicking  breast  phantoms, 
electric  field  measurements,  electrical-to-acoustical  transduction,  creating  a  3-D 
electromagnetic  stimulation  system,  and  validating  the  effectiveness  of  the 
multifrequency  adaptive  and  robust  technique  for  image  reconstruction. 


II.  Body 

II.  1  Experimental  Progress 

(a)  Initial  Phantom  Resonance  Experiments 

Experimentation  designed  to  verify  acoustic  resonance  of  the  tumor  phantom  by 
electromagnetic  excitation  is  underway.  The  experiment  setup,  shown  schematically  in 
Figure  1,  consists  of  an  electromagnetic  waveguide,  two  ultrasonic  transducers  (one  for 
transmission  ant  the  other  for  sensing),  and  a  Plexiglas  cylinder  that  contains  the  de¬ 
ionized  water  and  tumor  phantom.  The  cylinder  height  is  10.3  inches  and  the  inner 
diameter  is  7.25  inches.  The  transducer  used  in  the  experiment  is  a  0.5  inch  diameter 
piezoelectric  immersion  transducer  from  Olympus  NDT  (model  V303).  It  has  a  center 
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frequency  of  0.91  MHz  and  60.275%,  -6  dB  bandwidth.  At  1  MHz,  the  half  angle 
beam  width  is  approximately  3.5°.  Reciprocity  calibrations  are  being  performed  in  de¬ 
ionized  water.  The  waveguide  used  is  from  Penn  Engineering  Components  (model  WR 
187)  with  a  frequency  range  of  3.95-5.85  GHz.  In  the  ongoing  experiments,  a  1  cm 
diameter  sample  of  the  phantom  is  suspended  in  the  center  of  the  cylinder,  submerged  in 
de-ionized  water,  and  aligned  with  the  ultrasonic  transducers.  The  transmitting  transducer 
is  mounted  to  the  side  of  the  cylinder  with  the  transducer  face  exposed  to  the  de-ionized 
water.  The  receiving  transducer  is  similarly  aligned  and  its  axis  makes  a  180-degree  with 
respect  to  the  first  transducer.  The  electromagnetic  waveguide  provides  excitation  from 
the  bottom  of  the  cylinder.  Electric  field  probes  are  inserted  through  a  port  in  the  top  of 
the  cylinder  and  are  used  to  measure  total  field  intensity  along  the  cylinder  axis.  Two 
probes  are  needed:  one  to  measure  the  radial  field  intensity  and  another  to  determine  the 
axial  field  strength.  The  radial  field  probe  is  shown  in  Figure  4.  The  goal  of  the 
preliminary  setup  is  to  achieve  a  uniform  field  within  a  3  cm  radius  around  the  phantom. 
The  frequency  of  the  434  MHz  radio  frequency  signal  can  be  amplitude  modulated  over  a 
range  on  the  order  of  400  kHz  to  4  MHz,  which  is  the  expected  range  the  resonant 
frequency  of  the  phantom.  The  acoustic  signal  from  the  phantom  is  measured  by  the 
ultrasonic  transducer  and  is  recorded  by  the  NI  data  acquisition  system. 

Initially,  the  tumor’s  resonant  frequency  is  established  acoustically.  Referring  to 
Figure  1,  a  swept  sinusoidal  signal  generated  emanating  from  Port  1  of  an  RF  network 
analyzer  acts  as  a  source  which  drives  one  acoustic  transducer  while  the  second  acoustic 
transducer,  used  to  receive  the  acoustic  signal,  drives  port  2.With  the  tumor  absent  from 
the  water-filled  tank,  the  received  signal  consists  of  the  acoustic  frequency  characteristics 
of  the  tank  itself,  and  this  frequency  response  is  used  as  a  reference  to  calibrate  the 
network  analyzer.  The  same  sweep  is  then  performed  with  the  tumor  phantom  present 
and  the  response  now  contains  the  resonant  characteristics  of  the  phantom.  A  typical 
frequency  response  is  shown  in  Figure  2  where  it  is  seen  that  the  resonant  frequency 
occurs  between  around  1  MHz. 


Figure  1:  Test  setup  for  determining  the  resonant 
frequency  of  the  tumor  phantom. 
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Figure  2:  A  typical  frequency  response  of  acoustic  resonance. 

For  TAI  systems  the  acoustic  pressure  wave  is  generated  as  the  tumor  is  heated 
via  application  of  EM  energy.  Having  established  the  resonant  frequency  of  the  phantom 
using  the  acoustic  technique  described,  a  signal  with  this  frequency  amplitude  modulated 
at  434  MHz  EM  carrier  signal  is  then  used  to  drive  the  waveguide  located  at  the  bottom 
of  the  tank  via  an  RF  power  amplifier.  The  acoustic  transducer  is  then  used  to  sense  the 
acoustic  pressure  wave  generated  by  the  TAI  process. 

(b)  Phantom  Aging  Characteristics 

While  carrying  out  these  measurements,  it  became  apparent  the  large  number  of 
variables  involved  in  characterizing  the  resonance  of  the  phantom,  one  being  the  affect  of 
aging  of  the  phantom  over  time.  Subsequently,  some  time  was  focused  to  characterize 
the  phantom  aging.  One  of  the  aspects  explored  was  the  variation  of  the  resonant 
frequency  over  time.  Figure  3  shows  results  based  on  the  experimental  setup  described  in 
the  previous  section  and  shown  in  Figure  1.  While  the  sample  size  consists  of  only  3 
phantoms,  the  general  trend  seems  to  be  that  the  resonant  frequency  initially  increases 
before  it  decreases  over  time.  More  specimens  need  to  be  studied  to  confirm  this  trend. 
In  Figure  3,  phantom  c  exhibited  no  apparent  resonance  in  week  one.  This  was  attributed 
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to  phantom  placement  in  the  cylinder.  The  setup  was  consequently  altered  in  a  way  that 
would  insure  repeatable  phantom  placement.  Phantom  aging  is  important  to  quantify  in 
order  to  minimize  experimental  error  over  time.  Other  properties  are  currently  being 
examined  including  density  and  speed  of  sound. 
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Figure  3:  Resonant  frequency  of  phantom  over  time.  Note  no  resonance  was  found  for  phantom  c  in 
the  first  week. 

(c)  Electric  Field  Measurements 

Accurate  modeling  of  both  the  acoustical  and  electrical  behavior  of  the  system 
forms  an  essential  component  of  this  research  effort.  Of  special  concern  is  the  ability  to 
model  electrical-to-acoustical  transduction  aspect  of  this  system.  The  acoustic  pressure 
wave  to  be  measured  is  generated  when  the  modulated  electrical  signal  heats  the  phantom 
causing  it  to  mechanically  undulate.  This  heating  occurs  as  a  result  of  the  electric  field 
intensity  in  the  location  of  the  conducting  phantom.  As  a  consequence  not  only  must  the 
electrical  fields  present  in  this  system  be  accurately  modeled,  but  these  fields  must  also 
be  verified  so  that  the  proper  heating  parameters  may  be  used  as  input  to  the  acoustic 
model.  Initially  an  off-the-shelf,  low-profile,  high-resolution  probe  that  measures  the 
electric  potential  (Carsten  Associates  Model  E-601)  was  used.  Once  the  potential  is 
known  a  differencing  procedure  can  be  used  to  determine  the  field.  It  was  found,  however, 
that  this  method  proved  inadequate  for  the  present  system.  In  particular,  the  common 
mode  signal  measured  far  outweighed  the  difference  signal  and  an  accurate  estimate  of 
the  electric  field  could  not  be  determined.  Further,  the  metal  shaft  of  the  probe,  though 
small  in  diameter,  was  distorting  the  field  in  the  tank.  To  remedy  these,  two  custom 
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probes,  one  used  to  measure  the  axial  field,  another  used  to  measure  the  radial  field  were 
designed.  Since  the  new  probes  measure  the  electric  field  directly  no  differencing  is 
required.  Furthermore  the  new  probes  are  design  with  a  shaft  that  is  impedance-matched 
to  the  impedance  of  de-ionized  water  thus  making  the  probe  transparent  except  of  course 
(unavoidably)  for  the  immediate  region  being  measured.  To  date  the  radial  field  probe 
has  been  obtained  and  is  shown  in  Figure  4.  Shown  also  in  Figure  4  is  a  Transverse 
Electro-Magnetic  (TEM)  test  cell  used  for  probe  calibration.  One  end  of  the  TEM  cell  is 
excited  with  a  434  MHz  signal  while  the  other  end  is  terminated  in  a  50-ohm  load.  Since 
the  dimensions  of  the  cell  are  precisely  known,  the  electric  field  in  the  cell’s  center  is  also 
known.  This  center  portion  of  the  TEM  cell  is  filled  with  de-ionized  water  in  which  the 
probe  is  inserted,  and  the  resulting  voltage  is  measured  using  a  50-ohm  oscilloscope. 
Once  the  axial  filed  probe  is  completed,  an  accurate  field  map  of  the  tank  along  the  axis 
where  the  phantom  will  be  placed  will  be  obtained.  Once  again  this  data  will  be  used  to 
verify  the  electromagnetic  model  and  link  the  electrical  and  acoustic  models. 


Figure  4:  Radial  field  probe  and  probe  calibration  cell. 

(d)  Electrical-to- Acoustical  Transduction 

To  date  two  problems  have  been  encountered  with  regard  to  accurately  measuring 
the  RTAI  effect  and  are  currently  being  remedied.  The  first  regards  drift  of  the  acoustic- 
acoustic  resonance  measurements.  Since  the  acoustic  transducers  possess  narrow 
beamwidths,  accurate  placement  of  the  phantom  is  required  for  a  stable  resonance. 
Secondly  the  acoustic  resonance  varies  with  the  temperature  of  the  water.  Initially  the 
acoustic  and  electrical  portions  of  the  experimental  procedure  where  performed 
consecutively.  It  has  been  determined  that  monitoring  the  acoustic  resonance  while  the 
electrical  measurements  are  taken  is  essential.  Furthermore,  to  date  the  electrical  (RF) 
power  level  has  been  kept  low  with  a  total  drive  power  of  less  than  four  watts.  This 
requires  that  the  measurement  bandwidth  be  kept  small  so  as  to  not  allow  excessive  noise 
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in  the  measurements.  A  new  RF  power  amplifier  capable  of  producing  a  drive  power  of 
up  to  100- watts  has  been  obtained.  This,  along  with  accurate  placement  of  the  phantom’s 
location  will  allow  for  a  robust,  stable  measurement  of  the  RTAI  effect. 


II.  2  Electromagnetic  Stimulation 

Use  of  the  RTAI  concept  as  a  real  imaging  application  requires  extending  the 
approach  to  a  three-dimensional  system.  As  the  experimental  portion  of  present  phase  of 
the  research  effort  continues  so  does  the  modeling  for  the  next  phase  of  the  effort  which 
involves  three-dimensional  modeling,  simulation,  and  measurements. 


Figure  5  illustrates  conceptually  how  the  required  uniform  electric  field  for  a 
realistic  breast  phantom  can  be  achieved  in  a  simple  and  practical  manner.  Consider  the 
situation  shown  in  Figure  5(a)  which  illustrates  a  dielectric  sphere  of  arbitrary 
permittivity  in  free  space.  It  is  well  known  in  electromagnetic  theory  that  a  sinusoidally- 
distributed  surface  charge  distribution,  ps  -  po  cos  0 ,  with  po  as  a  constant,  placed  on  the 

surface  of  the  sphere  results  in  a  perfectly  uniform  electric  field  inside  the  dielectric 
sphere  as  shown.  By  invoking  the  theory  of  images  a  Perfect  Magnetic  Conductor  (PMC), 
which  represents  a  good  model  for  the  chest  wall,  can  be  placed  as  shown  without 
altering  the  electric  field  in  the  sphere.  In  this  way  a  uniform  electric  field  in  a 
hemispherical  region  that  models  the  human  breast  can  be  obtained.  Obtaining  such  a 
surface  charge  distribution  is,  in  general  not  possible,  though  it  can  be  approximated  in  a 
pointwise  fashion  as  illustrated  in  Figure  5(b).  There  it  is  seen  that  small  (subresonant) 
patches  can  be  placed  on  the  breast  so  as  to  approximate  a  uniform  electric  field  inside 
the  breast.  These  (metal)  patches  would  be  coated  with  a  non-conducting  material  and 
hence  modeled  as  capacitively  coupled  as  shown.  The  concept  is  akin  to  the  approach 
used  to  obtain  a  uniform,  circularly  polarized  transverse  magnetic  field  from  discrete  wire 
cage  structure  in  Magnetic  Resonance  Imaging  (MRI)  applications.  Though  the  perfectly 
uniform  field  would  be  obtained  only  for  the  case  of  an  electrostatic  charge  distribution, 
it  is  not  unreasonable  to  expect  a  quasi-static  solution  which  reasonably  approximates  the 
ideal,  static  situation  under  appropriate  conditions,  and  this  is  currently  being  explored. 
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Figure  5:  Motivation  for  achieving  a  uniform  field  within  a  hemispherical  region  using  small, 
capacitively  coupled  patches. 
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To  illustrate  the  basic  idea,  Figure  6  shows  a  pointwise  electrostatic  approximation  to  the 
continuous  case  using  31  point  charges  uniformly  distributed  over  the  breast  phantom 
surface  as  indicated  in  Figure  6(b).  The  magnitudes  of  the  point  charges  are  chosen  so 
that  they  equal  the  value  of  the  continuous  charge  distribution  at  that  point.  This  would 
serve,  for  example,  as  a  good  initial  guess  for  a  numerical  optimization  procedure  which 
would  give  the  optimum  charge  distribution  that  provide  the  most  uniform  electric  field 
magnitude  in  the  breast.  Figure  6(c)  shows  the  (unoptimized)  electrostatic  field  intensity 
over  several  planes  in  the  region  to  be  imaged.  It  is  seen  that  even  for  this  simple  case 
that  the  field  intensity  varies  by  less  than  1  dB  in  the  central  region  to  be  imaged.  This 
approach  is  currently  being  developed  analytically  as  well  as  being  modeled  using  two 
commercially  available  software  packages. 


Figure  6:  Simulated  uniform  electric  field  (quasi  TEM)  for  a  hemispherical  shell  excited  by  small 
patches. 


II.  3  Adaptive  Image  Formation  Algorithms 

One  of  the  key  challenges  encountered  by  thermal  acoustic  imaging  (TAI)  is  to 
develop  accurate  and  robust  image  reconstruction  methods.  Delay-and-sum  (DAS) 
beamforming  is  a  widely  used  reconstruction  algorithm  for  TAI.  However,  this  data 
independent  approach  tends  to  suffer  from  poor  resolution  and  high  sidelobe  level 
problems.  On  the  other  hand,  data  adaptive  approaches,  such  as  the  recently  introduced 
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robust  Capon  beamforming  (RCB)  [1]  method,  can  have  much  better  resolution  and 
much  better  interference  rejection  capabilities  than  their  data  independent  counterpart. 

We  propose  a  multifrequency  adaptive  and  robust  technique  (MART)  based  on 
RCB  for  TAI,  where  a  multiple  frequency  source  is  employed  instead  of  a  single 
frequency  source,  which  is  generally  used  by  other  research  groups.  MART  can  offer 
higher  signal-to-ratio  (SNR)  and  higher  imaging  contrast  than  its  single  frequency 
counterpart,  which  we  refer  to  as  the  single-frequency  adaptive  and  robust  technique 
(SART),  since  much  more  information  about  the  human  breast  can  be  harvested  from  the 
multiple  frequencies.  Furthermore,  the  interference  due  to  inhomogeneous  breast  tissue 
can  be  suppressed  more  effectively  since  more  information  about  the  breast  tissue  can  be 
used  by  the  RCB  algorithm.  MART  is  a  three-stage  time-domain  signal  processing 
algorithm.  In  Stage  I,  RCB  is  used  to  estimate  the  thermal  acoustic  responses  from  the 
focal  points  within  the  breast  for  each  stimulating  frequency.  Then  in  Stage  II,  a  scalar 
acoustic  waveform  at  each  focal  point  is  estimated  based  on  the  response  estimates  for  all 
frequencies  from  Stage  I.  Finally,  in  Stage  III,  the  positive  peak  and  the  negative  peak  of 
the  estimated  acoustic  waveform  at  each  grid  location  are  determined,  and  the  peak-to- 
peak  difference  is  computed  and  referred  to  as  the  image  intensity. 

To  validate  the  effectiveness  of  the  MART,  we  develop  a  2-D  inhomogeneous 
breast  model,  which  includes  skin,  breast  fatty  tissues,  glandular  tissues,  and  the  chest 
wall.  Small  tumors  are  set  below  the  skin.  The  finite-difference  time-domain  (FDTD) 
method  is  used  to  simulate  the  electromagnetic  field  inside  the  breast  tissues  [2,  3],  The 
specific  absorption  rate  (SAR)  distribution  is  calculated  based  on  the  simulated 
electromagnetic  field  [4,  5].  Then  FDTD  is  used  again  to  simulate  the  propagation  of  the 
thermal  acoustic  waves  [6,  7].  In  the  following  example,  the  thermal  acoustic  signals  are 
simulated  based  on  the  aforementioned  2-D  model  for  multiple  stimulating  frequencies 
from  200-800  MHz  with  frequency  step  100  MHz.  Two  small  1.5-mm-diameter  tumors 
are  set  inside  the  breast  model.  Their  locations  are  at  (X=70  mm,  Y=60  mm)  and  (X=75 
mm,  Y=62.5  mm).  The  distance  between  the  two  tumors  is  4  mm.  For  comparison 
purposes,  the  DAS  method  is  applied  to  the  same  data  set.  We  also  present  the  imaging 
results  obtained  by  SART  at  different  stimulating  frequencies.  Figure  7(a)  and  (b)  show 
the  imaging  results  for  MART  and  DAS,  respectively.  The  two  tumors  are  seen  clearly  in 
the  MART  image,  and  the  sizes  and  the  locations  of  the  two  tumors  are  accurate.  On  the 
other  hand,  the  DAS  image  contains  much  clutter  and  the  two  tumors  cannot  be  separated 
clearly.  Figure  7(c)  and  (d)  show  the  imaging  results  for  SART  at  the  stimulating 
frequencies  f= 200  and  800  MHz,  respectively.  The  tumors  can  be  seen  in  both  SART 
images,  but  more  clutters  show  up  in  Figure  7(c)  and  (d)  than  in  Figure  7(a),  and  the  sizes 
of  the  two  tumors  in  Figure  7(d)  are  larger  than  the  true  ones. 
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Figure  7:  Reconstructed  images  obtained  via  (a)  MART,  (b)  DAS,  (c)  SART  at  stimulating  frequency 
f=200  MHz,  (d)  SART  at  stimulating  frequency  f=800  MHz. 


III.  Key  Research  Accomplishments 

•  Phantom  development  was  continued  and  properties  as  a  function  of  aging  were 
explored 

•  Equipment  needed  for  RF  excitation  and  electric  field  measurements  was 
identified  and  purchased  (please  see  Table  1  for  details) 

•  Experimental  setup  to  detect  the  thermal  acoustic  waves  from  the  simulated  tumor 
was  formulated 

•  A  3-D  electromagnetic  stimulation  system  is  being  simulated  and  developed 

•  A  2-D  inhomogeneous  breast  model,  which  includes  skin,  breast  fatty  tissues, 
glandular  tissues,  and  the  chest  wall,  was  developed  to  validate  the  effectiveness 
of  the  MART 
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Vendor 

Model  number 

Equiment  purchased 

Qt. 

National  Instruments 

779079-02 

NI  PXI-5671  RF  Signal  Generator 

1 

Bruce  Carsten  Assoc. 

EFP  200R 

Calibrated  radial  probe 

1 

EFP  200A 

Calibrated  axial  probe 

1 

Mini-Circuits 

ZHL-100W-52-S 

Broadband  Amplifier 

1 

Acopian 

W530MT13 

Power  supply 

1 

Table  1:  Equipment  purchased  this  year. 


IV.  Reportable  Outcomes 

B.  Guo,  J.  Li,  H.  Zmuda,  and  M.  Sheplak,  “Multi-Frequency  Microwave  Induced 
Thermal  Acoustic  Imaging  for  Breast  Cancer  Detection,”  IEEE  Transactions  on 
Biomedical  Engineering,  Vol.  54,  No.  11,  pp.  2000-2010,  November  2007. 

Y.  Xie,  B.  Guo,  J.  Li,  G.  Ku,  and  L.V.  Wang,  “Adaptive  and  Robust  Techniques  (ART) 
for  Thermoacoustic  Tomography,”  to  appear  in  IEEE  Transactions  on  Biomedical 
Engineering. 

V.  Conclusion 

Preliminary  experiments  were  conducted  to  verify  the  feasibility  of  sensing  the 
acoustic  pressure  from  a  tissue  phantom  excited  at  resonance.  A  tissue  phantom  that 
mimics  the  electrical  and  mechanical  properties  of  a  tumor  for  use  in  the  excitation 
system  has  been  designed.  The  preliminary  phantom  development  focuses  on  matching 
dielectric  properties,  specifically  permittivity  and  conductivity.  The  phantom  was  made 
from  TX-151  powder  (from  Oil  Research  Center),  water,  cane  sugar  and  potassium 
chloride.  An  iterative  process  of  measuring  the  dielectric  constants  and  manipulating  the 
cane  sugar  and  potassium  chloride  concentration  was  conducted  until  a  suitable  phantom 
was  created. 

The  initial  experiment  has  been  designed  to  verify  acoustic  resonance  of  the 
tumor  phantom  by  electromagnetic  excitation.  The  experiment  consists  of  an 
electromagnetic  waveguide,  a  single  ultrasonic  transducer  for  sensing,  and  a  Plexiglas 
cylinder  that  contains  the  deionized  water  and  tumor  phantom. 

Currently  the  properties  of  the  tumor  phantom  (dielectric  constant,  conductivity, 
and  mechanical  resonant  frequency)  are  being  evaluated  as  a  function  of  temperature  and 
time  (age  of  the  phantom).  Preliminary  results  suggest  that  the  positioning  of  the  tumor 
phantom  within  the  Plexiglas  tank  strongly  influences  measured  results.  These 
characteristics  are  being  quantified  so  that  repeatable  results  can  be  obtained.  Since  this 
proof-of-concept  experimental  setup  uses  only  one  acoustic  transducer  and  one 
waveguide  (EM)  exciter,  such  sensitivity  is  not  unexpected.  The  next  generation  system 
which  uses  multiple  EM  exciters  and  multiple  acoustic  receivers  should  reduce  these 
sensitivities. 

We  have  devised  Adaptive  and  Robust  Methods  of  Reconstruction  (ARMOR)  for 
thermoacoustic  tomography.  ARMOR  is  robust  to  the  amplitude  and  phase  distortions  in 
the  recorded  signals  caused  by  the  acoustic  heterogeneity  of  biological  tissues.  Examples 
based  on  a  numerically  simulated  2-D  breast  model  and  two  sets  of  experimentally 
measured  data  from  human  mastectomy  specimens  demonstrate  the  excellent 
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performance  of  ARMOR:  high  resolution,  low-side  lobe  level,  and  much  improved 
interference  suppression  capability. 
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Abstract — Microwave-induced  thermal  acoustic  imaging  (TAI) 
is  a  promising  early  breast  cancer  detection  technique,  which  com¬ 
bines  the  advantages  of  microwave  stimulation  and  ultrasound 
imaging  and  offers  a  high  imaging  contrast,  as  well  as  high  spatial 
resolution  at  the  same  time.  A  new  multifrequency  microwave-in¬ 
duced  thermal  acoustic  imaging  scheme  for  early  breast  cancer 
detection  is  proposed  in  this  paper.  Significantly  more  information 
about  the  human  breast  can  be  gathered  using  multiple  frequency 
microwave  stimulation.  A  multifrequency  adaptive  and  robust 
technique  (MART)  is  presented  for  image  formation.  Due  to  its 
data-adaptive  nature,  MART  can  achieve  better  resolution  and 
better  interference  rejection  capability  than  its  data-independent 
counterparts,  such  as  the  delay-and-sum  method.  The  effective¬ 
ness  of  this  procedure  is  shown  by  several  numerical  examples 
based  on  2-D  breast  models.  The  finite-difference  time-domain 
method  is  used  to  simulate  the  electromagnetic  field  distribution, 
the  absorbed  microwave  energy  density,  and  the  thermal  acoustic 
field  in  the  breast  model. 

Index  Terms — Breast  cancer  detection,  finite-difference  time-do- 
main  (FDTD)  methods,  multifrequency  adaptive  and  robust 
technology  (MART),  robust  capon  beamforming  (RCB),  thermal 
acoustic  imaging  (TAI). 


I.  Introduction 

BREAST  cancer  is  the  most  common  nonskin  malignancy 
in  women  and  the  second  leading  cause  of  female  cancer 
mortality  [1].  There  are  over  200  000  new  cases  of  invasive 
breast  cancer  diagnosed  each  year  in  the  U.  S.,  and  one  out  of 
every  seven  women  in  the  U.S.  will  be  diagnosed  with  breast 
cancer  in  their  life  time  (the  American  Cancer  Society,  2006) 
and  early  diagnosis  is  key  to  surviving  breast  cancer  [2].  Mi¬ 
crowave  imaging  is  a  method  for  early  breast  cancer  detection 
[3]-[9],  which  exploits  the  significant  contrast  in  dielectric 
properties  between  normal  and  cancerous  tissues  [10]— [12]. 
However,  it  is  difficult  for  microwave  imaging  techniques  to 
achieve  good  (submillimeter)  spatial  resolution  because  of 
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Fig.  1.  Model  of  microwave-induced  TAI  for  breast  cancer  detection. 


its  long  wavelength  [13].  Ultrasound  is  another  option  which 
offers  a  high  spatial  resolution  because  of  its  short  acoustic 
wavelength  [14]— [16] .  However,  the  contrast  in  acoustic  prop¬ 
erties  between  normal  and  tumor  tissues  is  very  small  due  to 
both  being  soft  tissues. 

Microwave-induced  thermal  acoustic  imaging  (TAI)  com¬ 
bines  the  advantages  of  microwave  stimulation  and  ultrasound 
imaging  [13],  which  offers  a  high  imaging  contrast  (due  to  the 
significantly  different  dielectric  properties  of  tumor  and  normal 
breast  tissues),  as  well  as  high  spatial  resolution  (due  to  the  low 
propagation  velocity  or  the  short  wavelength  of  acoustic  waves 
in  biological  tissues)  at  the  same  time.  To  use  microwave-in¬ 
duced  TAI  techniques  for  breast  cancer  imaging,  a  microwave 
source  with  a  short  duration  time  is  used  to  irradiate  the  breast, 
as  shown  in  Fig.  1.  The  normal  breast  tissues,  as  well  as  tumors, 
absorb  microwave  energy  and  emanate  thermal  acoustic  waves 
by  thermoelastic  expansion.  It  is  well-known  that  malignant 
breast  tissue  has  a  higher  water  content  [1],  [10],  [12],  [17], 
with  a  much  higher  conductivity  than  normal  breast  tissues 
(with  low  water  content).  As  a  result,  the  microwave  energy  ab¬ 
sorbed  by  tumor  and  normal  breast  tissues  will  be  significantly 
different  and  a  stronger  acoustic  wave  will  be  produced  by  the 
tumor.  The  acoustic  waves  generated  in  this  manner  carry  the 
information  about  the  microwave  energy  absorption  properties 
of  the  tissues  under  irradiation.  The  thermal  acoustic  waves 
propagate  out  of  the  breast  and  are  recorded  by  an  acoustic 
sensor  array  placed  around  the  breast.  The  tumor  locations  can 
be  accurately  determined  since  the  received  acoustic  signals 
from  the  malignant  tumors  are  at  higher  levels,  hence  aiding 
image  construction. 

During  the  last  decade,  several  research  groups  have  been 
working  on  the  microwave-induced  TAI  of  biological  tissues 
[1 8]— [23] .  The  microwave  frequency  used  ranges  from  400  MHz 
[22]  to  3  GHz  [13].  Image  reconstruction  algorithms  include  the 
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widely  used  delay-and-sum  (DAS)  method  [20],  [23],  the  fre¬ 
quency-domain  inverse  method  [24],  [25],  and  the  time-domain 
inverse  method  [13],  [19]. 

Microwave-induced  TAI  does,  however,  present  several 
challenges.  First,  the  human  breast  is  large  in  size,  usually 
has  an  irregular  shape  if  not  compressed,  and  is  covered  with 
a  2-mm-thick  skin  with  dielectric  properties  significantly 
different  from  the  normal  breast  tissues.  Moreover  the  breast 
tissue  is  far  from  homogeneous  because  the  dielectric  proper¬ 
ties  of  glandular  tissue  are  different  from  that  of  breast  fatty 
tissue.  All  these  factors  make  it  difficult  to  approximate  the 
back  propagation  properties  of  thermal  acoustic  signals  inside 
the  breast.  Due  to  the  slow  acoustic  wave  propagation  speed  or 
short  wavelength  in  biological  tissues,  the  errors  on  the  order 
of  millimeters  in  determining  the  acoustic  signal  propagation 
path  lengths  will  severely  degrade  the  image  quality. 

In  this  paper,  a  multifrequency  microwave-induced  TAI 
system  is  proposed  which  remedy  the  problems  mentioned 
above.  Instead  of  using  a  single  frequency  microwave  source, 
as  generally  done  by  other  research  groups  in  this  field,  here 
a  multiple  frequency  source  is  used,  since  the  desired  thermal 
acoustic  signals  can  be  induced  by  microwave  sources  operating 
at  a  wide  range  of  frequencies.  We  show  in  this  paper  that  the 
rich  information  collected  from  the  multifrequency  stimulation 
can  help  mitigate  the  challenges  mentioned.  The  multifre¬ 
quency  microwave-induced  thermal  acoustic  signals  will  offer 
higher  signal-to-noise  ratio  (SNR)  and  higher  imaging  contrast 
than  single-frequency  microwave-induced  thermal  acoustic  sig¬ 
nals  since  much  more  information  about  the  human  breast  can 
be  harvested  from  the  multiple  stimulating  frequencies  within 
the  microwave  frequency  band.  Furthermore,  the  interference 
due  to  inhomogeneous  breast  tissue  can  be  suppressed  more 
effectively  when  multifrequency  microwave-induced  thermal 
acoustic  signals  are  used  for  image  reconstruction  since  more 
information  about  breast  tissue  can  be  used  by  the  adaptive 
image  reconstruction  algorithms. 

Another  challenge  encountered  by  microwave-induced  TAI 
is  the  need  to  develop  accurate  and  robust  image  reconstruction 
methods.  DAS  is  a  widely  used  reconstruction  algorithm  in 
medical  imaging.  This  method  is  data-independent  and  tends  to 
suffer  from  poor  resolution  and  high  sidelobe  level  problems. 
Data-adaptive  approaches,  such  as  the  recently  introduced 
robust  Capon  beamforming  (RCB)  [26],  [27]  method,  can  have 
much  better  resolution  and  much  better  interference  rejection 
capabilities  than  their  data-independent  counterpart.  Several 
medical  imaging  algorithms  [4],  [5],  [28],  [29]  based  on  RCB 
have  been  developed  and  used  for  microwave  imaging  and 
thermal  acoustic  imaging.  Good  performances  of  these  algo¬ 
rithms  have  been  reported. 

We  present  a  multifrequency  adaptive  and  robust  technique 
(MART)  based  on  RCB  for  multifrequency  microwave-induced 
TAI.  There  are  three  stages  in  our  MART.  In  Stage  I,  RCB  is 
used  to  estimate  the  thermal  acoustic  responses  from  the  grid 
locations  within  the  breast  for  each  stimulating  microwave 
frequency.  Then,  in  Stage  II,  a  scalar  acoustic  waveform  at 
each  grid  location  is  estimated  based  on  the  response  estimates 
for  all  stimulating  frequencies  from  Stage  I.  Finally,  in  Stage 
III,  the  positive  peak  and  the  negative  peak  of  the  estimated 


acoustic  waveform  at  each  grid  location  are  determined,  and 
the  peak-to-peak  difference  is  computed  and  referred  to  as  the 
image  intensity. 

To  validate  the  effectiveness  of  the  proposed  algorithm,  we 
develop  a  2-D  inhomogeneous  breast  model,  which  includes 
skin,  breast  fatty  tissues,  glandular  tissues,  and  the  chest  wall. 
Small  tumors  are  set  below  the  skin.  The  finite-difference  time- 
domain  (FDTD)  method  is  used  to  simulate  the  electromagnetic 
field  inside  the  breast  tissues  [30],  [31].  The  specific  absorp¬ 
tion  rate  (SAR)  distribution  is  calculated  based  on  the  simu¬ 
lated  electromagnetic  field  [32],  [33].  Then  FDTD  is  used  again 
to  simulate  the  propagation  of  the  microwave-induced  thermal 
acoustic  waves  [34],  [35]. 

The  remainder  of  this  paper  is  organized  as  follows.  In 
Section  II,  the  microwave  frequency  properties  of  human 
breast  are  described.  A  proper  microwave  frequency  band  for 
multifrequency  microwave-induced  TAI  is  also  given  in  this 
section.  MART  is  proposed  for  image  formation  in  Section  III. 
In  Section  IV,  2-D  electromagnetic  and  acoustic  breast  models 
are  developed.  The  electromagnetic  and  acoustic  simulation 
methods  are  also  presented  in  this  section.  Imaging  results 
based  on  numerical  examples  are  provided  in  Section  V. 
Section  VI  concludes  the  paper. 

II.  Microwave  Properties  of  Human  Breast 

A.  Cutoff  Frequency  of  Human  Breast 

In  a  microwave-induced  TAI  system,  the  biological  tissues 
should  be  heated  by  microwave  sources  in  a  uniform  manner, 
otherwise  thermal  acoustic  signals  will  be  induced  by  a  nonuni¬ 
form  microwave  energy  distribution,  resulting  in  images  diffi¬ 
cult  to  interrupt.  It  is  well-known  that  high-order  electromag¬ 
netic  field  modes  will  be  excited  in  a  media  if  the  microwave 
works  at  a  frequency  higher  than  a  cutoff  frequency  of  the  media 
[36],  and  the  microwave  energy  distribution  is  nonuniform  at 
high-order  modes  [37].  To  minimize  the  nonuniform  microwave 
energy  distribution  inside  the  breast  caused  by  the  high-order 
electromagnetic  modes,  the  microwave  source  should  work  at  a 
frequency  below  a  certain  cutoff  frequency. 

To  estimate  the  cutoff  frequency  for  the  human  breast,  we 
consider  the  simplified  breast  model  shown  in  Fig.  2(a)  con¬ 
sisting  of  a  semicircular  dielectric  waveguide  with  a  perfect 
magnetic  conductor  (PMC)  at  the  bottom  of  the  semicircle.  Re¬ 
call  that  the  tangential  components  of  the  magnetic  field  are  zero 
on  the  surface  at  the  PMC.  The  PMC  assumption  is  reasonable 
because  the  permittivity  of  the  chest  wall  (er  =  50)  is  much 
greater  than  that  of  the  normal  breast  tissues  (er  =  9).  In  cir¬ 
cular  dielectric  waveguide,  if  an  electromagnetic  mode  has  a 
field  distribution  whose  tangential  magnetic  field  components 
are  zero  at  the  center  line  of  the  circular  waveguide,  as  shown 
in  Fig.  2(b),  the  introduction  of  a  PMC  at  the  center  line  of  the 
circular  waveguide  will  not  significantly  change  the  boundary 
conditions  and,  hence,  will  not  significantly  alter  the  mode  dis¬ 
tribution.  The  modes  in  the  semicircular  dielectric  waveguide 
can,  thus,  be  estimated  by  determining  the  modes  in  a  corre¬ 
sponding  circular  dielectric  waveguide. 
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Perfect  Magnetic  Conductor  (PMC) 
(a) 


HE2 1  TM01 

Electric  field  -  Magnetic  field  - - 

(b) 

Fig.  2.  Simplified  breast  model,  (a)  Semicircular  dielectric  waveguide  with 
PEC,  and  (b)  corresponding  circular  dielectric  waveguide. 


The  dominant  mode  of  a  circular  dielectric  waveguide  is 
the  HE  11  mode,  the  cutoff  frequency  of  which  is  zero.  The 
electromagnetic  field  distribution  is  near  uniform  at  this  mode. 
The  dominant  mode  is  followed  by  the  TE01,  TM01,  and 
HE21  modes.  These  modes  are  degenerate,  and  have  a  cutoff 
frequency  given  by  [36] 


fc  = 


XoiCq 
2'Ka\Jer  —  1 


(1) 


where  Co  is  the  speed  of  light  in  free  space,  xoi  —  2.405  is  the 
first  root  of  the  Bessel  function  of  the  first  kind  of  order  zero 
(^o(Xoi)  =  0),  a  and  er  are  the  radius  and  average  permit¬ 
tivity  of  the  circular  dielectric  waveguide,  respectively.  TM01 
and  HE21,  as  well  as  the  interference  between  them,  satisfy  the 
zero  tangential  magnetic  field  component  condition  at  the  center 
line  of  the  circular  waveguide.  These  modes  can  also  exist  in  the 
semicircular  dielectric  waveguide.  By  substituting  the  parame¬ 
ters  of  the  breast  model  into  (1),  we  obtain  the  cutoff  frequency 
of  the  semicircular  breast  model  to  be 


fc  = 


2.405  •  3  x  108 
2t r  •  0.05  • 


=  812  MHz 


(2) 


where  we  have  used  a  —  5  cm  and  er  —  9  as  typical  values 
for  human  breast.  Consequently,  the  stimulating  microwave  fre¬ 
quency  for  the  TAI  system  should  be  below  812  MHz. 


B.  Microwave  Energy  Absorption  Properties  of  Breast  Tissues 
and  Tumor 

It  is  well-known  that  the  complex  relative  dielectric  properties 
of  a  medium  can  be  expressed  as 


<-V  & r 


TABLE  I 

Cole-Cole  Parameters  for  Biological  Tissues 


Tissue 

Breast 

Skin 

Muscle 

Tumor 

£oo 

2.5 

4.0 

4.0 

4.0 

a 

0.01 

0.0002 

0.2 

0.7 

Asx 

3.0 

32.0 

50.0 

50.0 

n  (ps) 

17.68 

7.23 

7.23 

7.0 

ChL\ 

0.1 

0.0 

0.1 

0.0 

As2 

15 

1100 

7000 

0 

T2  (ns) 

63.66 

32.48 

353.68 

N/A 

OL2 

0.1 

0.2 

0.1 

N/A 

As3 

5.0E4 

0 

1.2E6 

0 

T3  (MS) 

454.7 

N/A 

318.31 

N/A 

G3 

0.1 

n/a 

0.1 

N/A 

A^4 

2.0E7 

n/a 

2.5E7 

0 

T4  (ms) 

13.26 

n/a 

2.274 

N/A 

CK4 

0.0 

n/a 

0.0 

N/A 

where  e'r  is  the  relative  permittivity  and  s"  is  the  out-of-phase 
loss  factor  which  can  be  written  as 


with  a  being  the  total  conductivity,  £o  the  free  space  permit¬ 
tivity,  and  u  the  electromagnetic  frequency.  The  tissue  absorp¬ 
tion  property  of  the  electromagnetic  wave  energy  is 

Q{v)=1-cj\E{v)\2  (5) 

which  is  a  function  of  the  total  conductivity  and  the  electric  field 
inside  the  tissue.  If  we  assume  that  the  microwave  energy  distri¬ 
bution  is  uniform  inside  the  breast  in  a  TAI  system,  the  absorp¬ 
tion  of  the  microwave  energy  by  the  breast  is  characterized  by 
the  total  conductivity  of  the  breast  tissues 


ct(uj)  =  e"eo  u. 


(6) 


Hence,  instead  of  using  the  attenuation  coefficient  a,  as  used  in 
[23],  in  this  paper,  we  study  the  absorption  properties  of  breast 
tissues  using  the  total  conductivity  a. 

The  dielectric  properties  of  biological  tissues  can  be  accu¬ 
rately  modeled  by  the  Cole-Cole  equation  [38] 


K 

£r{ij)  =  Sqo  + 

i= 1 


A£j 

1  +  (jWTi)1-** 


+ 


^0 
jus  0 


(7) 


where  K  is  the  order  of  the  Cole-Cole  model,  is  the 
high-frequency  permittivity,  rv;  is  the  relaxation  time,  A Si  is 
the  pole  amplitude,  ai  (0  <  a.i  <  1)  is  a  measure  of  the  broad¬ 
ening  of  dispersion,  and  ao  is  the  static  ionic  conductivity.  The 
Cole-Cole  parameters  for  skin,  breast  fatty  tissue,  chest  wall 
(mainly  consisting  of  muscle),  as  well  as  tumor  are  listed  in 
Table  I  [39],  [40].  Because  we  cannot  find  the  values  specific  to 
the  tumor,  the  dielectric  properties  of  the  tumor  is  approximated 
using  a  Debye  model  [3],  [41],  which  is  a  special  case  of  the 
Cole-Cole  model. 

Substituting  (7)  into  (6),  we  obtain  the  total  conductivity  of 
the  breast  tissue  as  follows: 


K 


r(cj)  =  -imag  (  £oo  +  E,  i 


A  £i 


Cl1  + 


+ 


£0CU 

(8) 


(3) 
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Fig.  3.  Total  conductivity  of  normal  breast  tissues  and  tumor  as  a  function  of 
frequency. 


Fig.  4.  Ratio  of  conductivity  between  tumor  and  normal  breast  tissue  as  a  func¬ 
tion  of  frequency. 


which  is  a  function  of  the  stimulating  microwave  frequency, 
where  imag(-)  denotes  the  imaginary  part  of  the  complex 
relative  permittivity.  Fig.  3  gives  the  total  conductivity  of 
breast  fatty  tissue  and  tumor  over  a  frequency  band  from 
100-1000  MHz.  Note  that  the  total  conductivity  increases  with 
the  microwave  frequency,  which  means  that  more  microwave 
energy  is  absorbed  and  converted  to  heat  by  tissues  at  higher 
microwave  frequency  region,  or  in  other  words,  the  SNR  is 
higher  in  the  received  thermal  acoustic  signals  at  higher  stim¬ 
ulating  microwave  frequency  region.  On  the  other  hand,  the 
penetration  at  higher  microwave  frequencies  is  smaller  because 
the  tissues  are  lossy.  We  define  the  conductivity  ratio  between 
the  tumor  and  the  normal  breast  tissue  as 


r<r(w)  = 


^tumor(^) 
<T  breast  (^) 


(9) 


and  plot  it  as  a  function  of  frequency  in  Fig.  4.  A  high  con¬ 
ductivity  ratio  means  that  more  microwave  energy  is  absorbed 
and  converted  to  heat  by  tumor  than  by  normal  breast  tissues. 
In  other  words,  the  higher  the  conductivity  ratio,  the  higher 
the  imaging  contrast.  Fig.  4  shows  that  the  imaging  contrast 
is  higher  at  the  lower  microwave  frequency  region  because  the 
conductivity  ratio  decreases  with  the  frequency. 

These  microwave  energy  absorption  properties  of  breast 
tissues  and  tumor  motivate  us  to  consider  inducing  thermal 
acoustic  signals  with  different  microwave  frequencies.  By 


taking  into  account  the  aforementioned  cutoff  frequency  given 
in  (2),  we  choose  a  frequency  range  from  200-800  MHz.  The 
frequency  step  is  100  MHz,  with  a  total  of  seven  frequencies. 
Note  that  wideband  antenna  techniques  should  be  used  for 
the  practical  implementation  because  the  frequency  range  is 
wide.  However,  since  the  exciting  microwave  frequency  is 
stepped,  an  antenna  with  a  broad  instantaneous  bandwidth 
is  not  required.  Another  advantage  of  using  multiple  fre¬ 
quencies  for  stimulation  is  that  more  information  about  the 
inhomogeneous  breast  tissues  will  be  harvested  from  the  mul¬ 
tifrequency  microwave-induced  thermal  acoustic  signals.  The 
microwave  energy  distribution  inside  the  breast  model  is  not 
uniform  because  the  human  breast  is  inhomogeneous  media, 
and  thermal  acoustic  signals  will  be  induced  by  the  inhomo¬ 
geneous  energy  distribution.  These  thermal  acoustic  signals 
will  appear  as  clutter  in  the  resulting  images.  However,  the 
inhomogeneous  microwave  energy  distributions  are  different  at 
different  stimulating  frequencies  because  of  the  different  mi¬ 
crowave  wavelengths  in  breast  tissues.  When  a  multifrequency 
microwave  source  is  used  for  TAI,  the  thermal  acoustic  clutter 
induced  by  the  inhomogeneous  breast  tissues  can  be  suppressed 
by  our  adaptive  and  robust  imaging  algorithm. 

III.  Multifrequency  Adaptive  and  Robust  Technology 
(MART)  for  Breast  Cancer  Imaging 

We  consider  a  multifrequency  microwave-induced  TAI 
system  as  shown  in  Fig.  1,  where  an  acoustic  sensor  array  is 
arranged  on  a  semicircle  relatively  close  to  the  breast  skin.  The 
location  of  each  acoustic  sensor  is  rj  (j  =  1,  •  •  • ,  TV),  where 
N  is  the  number  of  acoustic  sensors.  Assume  that  M  —  7  mi¬ 
crowave  sources  with  different  frequencies  are  used  to  irradiate 
the  breast  model.  Let  Pij(t )  ( i  =  1,  •  •  • ,  M\  j  =  1 ,  •  •  • ,  TV; 
t  =  0 ,  •  •  • ,  T  —  1)  denote  the  thermal  acoustic  signal  induced 
by  the  ith  frequency  and  received  by  the  jth  acoustic  sensor, 
where  T  is  the  recording  time  which  is  sufficiently  long  to 
allow  acoustic  sensors  to  record  all  responses  from  the  breast. 
Our  goal  is  to  detect  the  tumor  by  reconstructing  an  image  of 
the  thermal  acoustic  response  intensity  J(r)  as  a  function  of 
scan  location  r  within  the  breast. 

A.  Data  Preprocessing 

Because  the  breast  skin,  breast  tissues,  chest  wall,  and  tumor 
absorb  the  microwave  energy  and  convert  the  energy  to  heat,  all 
of  them  produce  thermal  acoustic  signals.  The  received  thermal 
acoustic  waveforms  include  the  responses  from  the  tumor,  as 
well  as  from  other  healthy  breast  tissues.  The  thermal  acoustic 
signals  generated  by  the  skin  are  much  stronger  than  those  by  a 
small  tumor  because  of  the  high  conductivity  of  the  skin  and  the 
acoustic  sensors  being  very  close  to  the  skin.  We  must  remove 
the  skin  responses  to  enhance  the  tumor  responses.  Because  the 
distances  between  the  acoustic  sensors  and  the  nearest  breast 
skin  are  similar  to  one  another,  the  signals  recorded  by  var¬ 
ious  sensors  have  similar  skin  responses.  Hence,  we  can  remove 
the  skin  response  by  subtracting  out  a  fixed  calibration  signal 
from  all  received  signals.  This  calibration  signal  can  be  obtained 
simply  by  averaging  the  recorded  signals  from  all  channels. 

Let  Xij(t)  denote  the  signals  after  subtracting  out  the  calibra¬ 
tion  signal.  To  process  the  signals  coherently  for  a  focal  point 
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at  r,  we  align  the  signals  Xij(t)  by  time  shifting  each  signal  a 
number  of  samples  77  ( r).  The  discrete  time  delay  between  r 
and  the  jth  acoustic  sensor  can  be  calculated  as 


nA r)  = 


(10) 


where  |_7j  stands  for  rounding  to  the  greatest  integer  less  than 
7,  1 1 17  —  r||  is  the  distance  between  17  and  r,  c  is  the  velocity 
of  the  acoustic  wave  propagating  in  breast  tissues,  and  At  is  the 
sampling  interval,  which  is  assumed  to  be  sufficiently  small.  The 
time- shifted  signals  are  denoted  as 


r)  =  Xij  (t  +  nj( r)) ,  t  =  -nj( r),  ■  ■  ■  ,T  -  nj( r). 

(11) 

After  time  shifting,  the  acoustic  signals  from  the  imaging  lo¬ 
cation  r  are  aligned  so  that  they  all  start  approximately  from 
time  t  =  0  for  all  channels.  Now  the  aligned  signals  are  win¬ 
dowed  by 


Window(Z)  = 


1,  0  <  /  <  L  —  1 
0,  otherwise 


(12) 


to  isolate  the  signals  from  the  focal  point  at  r.  The  windowed 
signals  are  denoted  as  X{j(l,  r),  /  =  0,  •  •  • ,  L  —  1,  where  LAt  is 
the  approximate  duration  of  the  thermal  acoustic  pulse,  which 
can  be  determined  from  the  pulse  duration  of  the  pulsed  mi¬ 
crowave  source. 

Attenuation  exists  when  acoustic  waves  propagate  within  the 
breast.  This  attenuation  has  two  parts:  the  attenuation  due  to  the 
lossy  media  and  the  propagation  attenuation.  Thus,  the  attenu¬ 
ation  of  the  tumor  responses  at  various  channels  are  different 
because  of  the  different  distances  between  the  imaging  position 
r  and  the  acoustic  sensors.  For  the  2-D  case  considered  here,  the 
compensation  factor  for  the  jth  channel  is  given  by 


Fig.  5.  Data  cube  model.  In  Stage  I,  MART  slices  the  data  cube  for  each  fre¬ 
quency  index.  RCB  is  applied  to  each  data  slice  to  estimate  the  corresponding 
waveform. 


MART  is  a  three-stage  time-domain  signal  processing  algo¬ 
rithm.  In  Stage  I,  MART  slices  the  data  cube  corresponding  to 
each  frequency  index,  and  processes  each  data  slice  by  the  RCB 
to  obtain  the  thermal  acoustic  waveform  estimate  for  each  stim¬ 
ulating  frequency.  Then,  in  Stage  II,  a  scalar  waveform  is  es¬ 
timated  from  all  frequencies  based  on  the  waveform  estimates 
from  Stage  I.  Finally,  the  positive  peak  and  the  negative  peak 
of  the  estimated  thermal  acoustic  waveform  from  Stage  II  are 
found  in  Stage  III.  The  peak-to-peak  difference  is  calculated  as 
the  image  intensity  at  the  focal  point  at  r.  The  details  of  all  three 
stages  are  given  below. 

1 )  Stage  I:  In  Stage  I,  MART  approximates  the  data  model 
as 


yi(Z)  =  a iSi{l)  +  e;(Z)  (16) 


Kj(r)  =  exp  (a\\rj  -  r||)  •  ||rj  -  r||1/2  (13) 

where  the  first  term  of  the  right  side  of  (13)  compensates  for  the 
attenuation  due  to  the  lossy  media,  and  the  second  term  com¬ 
pensates  for  the  geometric  attenuation.  The  compensated  signal 
can  be  calculated  as 

Vi, r)  =  KAr)  -  %i(f  r|s  l  =  0,  -  •  ■  ,L  —  1.  (14) 


B.  Multifrequency  Adaptive  and  Robust  Technology  (MART) 

Without  loss  of  generality,  we  consider  imaging  at  a  generic 
location  r  only.  For  notational  convenience,  we  drop  the  depen¬ 
dence  of  yij(l,  r)  on  r,  and  simply  denote  it  as  Now  we 

consider  the  data  model 

Vi.j  (0  —  $i,j  (0  “I”  tii.  j  (/)  (15) 

where  represents  the  tumor  response  and  eij  (/)  rep¬ 

resents  the  residual  term,  which  includes  the  noise  and  inter¬ 
ference  from  breast  skin,  chest  wall,  and  other  responses.  The 
structure  of  the  data  model  is  a  data  cube  as  shown  in  Fig.  5. 


where  y;(Z)  =  ■  ■  ■ ,  2/i,jv(Z)]T  and  e;(Z)  = 

[ei,i(Z)>  •  •  •  j ei,iv(Z)]T-  The  scalar  waveform  .s.,(Z)  denotes 
the  thermal  acoustic  signal  generated  at  the  focal  location  r 
corresponding  to  the  7th  stimulating  frequency.  The  vector  a i  is 
referred  to  as  the  array  steering  vector,  which  is  approximately 
equal  to  l.vx  i  =  [1,  •  •  • ,  1]T  since  all  the  signals  have  been 
aligned  temporally  and  their  attenuation  compensated  for  in  the 
preprocessing  step.  The  residual  e^(Z)  is  the  noise  and  interfer¬ 
ence  term,  which  is  assumed  uncorrelated  with  the  signal. 

There  are  two  assumptions  made  to  write  the  model  given  in 
(16).  First,  the  steering  vector  is  assumed  to  vary  with  the  mi¬ 
crowave  frequency  (i)  but  nearly  constant  with  the  time  sample 
(/).  Second,  we  assume  that  the  thermal  acoustic  signal  wave¬ 
form  depends  only  on  the  microwave  frequency  (i)  but  not  on 
the  acoustic  sensor  (j ) .  The  truth,  however,  is  that  the  steering 
vector  is  not  exactly  known  as  it  changes  slightly  with  both  the 
stimulating  frequency  and  time  due  to  array  calibration  errors 
and  other  factors.  The  signal  waveform  can  also  vary  slightly 
with  both  the  stimulating  frequency  and  acoustic  sensor,  due  to 
the  inhomogeneous  and  frequency-dependent  medium  within 
the  breast.  The  two  aforementioned  assumptions  simplify  the 
problem  slightly.  They  cause  little  performance  degradations 
when  used  with  our  adaptive  and  robust  algorithm. 
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In  practice,  the  true  steering  vector  in  (16)  is  not  ljvxi-  We 
assume  that  the  true  steering  vector  a i  lies  in  the  vicinity  of  the 
assumed  steering  vector  a  =  1  jyx.i ,  and  the  only  knowledge  we 
have  about  a i  is  that 


|a;-a||2  <*,  (17) 

where  e\  is  a  user  parameter,  which  may  be  determined  de¬ 
pending  on  the  various  errors  discussed  previously. 

The  true  steering  vector  a^  can  be  estimated  via  the  following 
covariance  fitting  approach  of  RCB  [26],  [27] 


ma xSf  subject  to  Ry,  —  >  0 

%  5  1 

||a*  -  a||2  <  6!  (18) 

where  Sf  is  the  power  of  the  signal  s*(Z)  and 

Ry  =  t  y^y,:(Z)yf(/)  (19) 

J  1= 0 

is  the  sample  covariance  matrix.  The  above  optimization 
problem  can  be  solved  as  described  in  [26],  and  the  estimated 
true  steering  vector  is  denoted  here  as  a,;. 

To  obtain  the  signal  waveform  estimate,  we  pass  the  received 
signals  through  a  Capon  beamformer  [27],  [42].  The  weight 
vector  of  the  beamformer  is  determined  by  using  the  estimated 
steering  vector  a,  in  the  following  expression: 


r;1! 


Wi  = 


am-1 


•Yi 


(20) 


Then  the  estimated  signal  waveform  corresponding  to  the  ith 
stimulating  frequency  is 


k(l)  =  w/y;(7).  (21) 


By  repeating  the  aforementioned  process  for  i  =  1  through 
i  =  M,  we  obtain  the  complete  set  of  M  waveform  estimates 

s(Z)  =  [si(Z),  •  •  • ,  sm(1)]T  ■  (22) 


2)  Stage  II:  Since  the  stimulating  microwave  sources  with 
various  frequencies  are  assumed  to  have  the  same  power,  we  as¬ 
sume  that  the  thermal  acoustic  responses  from  the  tumor  at  dif¬ 
ferent  stimulating  frequencies  have  nearly  identical  waveforms. 
Note  that  the  thermal  acoustic  responses  induced  by  the  inho¬ 
mogeneous  microwave  energy  distribution  (due  to  the  inhomo¬ 
geneous  breast  tissues)  are  different  at  different  stimulating  fre¬ 
quencies.  This  means  that  the  elements  of  the  vector  s (Z)  are  all 
approximately  equal  to  an  unknown  scalar  signal  s(Z),  and  the 
noise  and  interference  term  can  be  assumed  uncorrelated  with 
this  signal.  In  Stage  II  of  MART,  we  adopt  the  data  model 


where  as  is  approximately  equal  to  Imxi-  However,  the 
“steering  vector”  may  again  be  imprecise,  and,  hence,  RCB  is 
needed  again. 

As  we  did  in  Stage  I,  we  assume  that  the  only  knowledge 
about  as  is  that 


||as  -as||2  <  e2,  (24) 

where  as  =  1m  xi  is  the  assumed  steering  vector,  and  62  is  a 
user  parameter.  Again,  the  true  steering  vector  as  can  be  esti¬ 
mated  via  the  covariance  fitting  approach 

max#2  subject  to  Rs  —  £2asaJ  >  0 

62  ,as 

||as  -  a, ||2  <  e2  (25) 


where  82  is  the  power  of  the  signal  s(l),  and 

r  »  =  2  j;s(osT(o 

'  1=0 


(26) 


is  the  sample  covariance  matrix. 

After  obtaining  the  estimated  steering  vector  as,  we  obtain 
the  adaptive  weight  vector  and  the  estimated  signal  waveform, 
respectively,  as 


Rs  1a8 

ajRUas 


(27) 


and 


s(t )  =  wTs  (t).  (28) 

3)  Stage  III :  Because  the  thermal  acoustic  pulse  is  usually 
bipolar:  a  positive  peak,  corresponding  to  the  compression 
pulse,  and  a  negative  peak,  corresponding  to  the  rarefaction 
pulse  [43],  we  use  the  peak-to-peak  difference  as  the  response 
intensity  for  the  imaging  location  r  in  the  third  stage  of  MART. 
Compared  with  other  energy  or  amplitude  based  response 
intensity  estimation  methods,  peak-to-peak  difference  can  be 
used  to  improve  imaging  quality  with  little  additional  compu¬ 
tation  costs. 

The  positive  and  negative  peak  values  of  the  estimated  wave¬ 
form  for  the  focal  location  r  will  be  searched  based  on  the 
estimated  waveform  (28)  obtained  in  Stage  II.  Because  of  the 
nonuniform  sound  speed  in  biological  tissues,  the  arrival  time 
of  the  acoustic  pulse  generated  at  location  r  cannot  be  calcu¬ 
lated  accurately.  However,  it  was  reported  in  [18]  that,  when 
the  heterogeneity  is  weak,  such  as  in  breast  tissues,  amplitude 
distortion  caused  by  multipath  is  not  severe.  We  assume  that  the 
original  peak  remains  a  peak  in  the  estimated  waveform,  and  the 
positive  and  negative  peak  values  of  the  thermal  acoustic  pulse 
can  be  searched  as 

P+  =  max  <  max  |s(Z)}  ,  0  1  (29) 

( /G[A  i .  A-2,]  J 


s(Z)  =  a  ss(Z)  +es(Z) 


(23) 


and 


P  =  min  <  min  {£(/)},  0 

(^G[Ai,A2] 


(30) 
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Fig.  6.  Breast  model  for  thermal  acoustic  simulation,  (a)  Model  for  electro¬ 
magnetic  simulation  and  (b)  model  for  acoustic  simulation. 


where  [Ai,  A2]  E  [0,1/]  is  the  searching  range.  Here  Ai  and 
A2  are  user  parameters,  and  the  details  on  how  to  choose  them 
can  be  found  in  [29]. 

After  the  positive  and  negative  peak  values  are  found,  the 
response  intensity  for  the  focal  point  at  location  r  is  given  as 

J(r)  =  P+  —  P~ .  (31) 


IV.  Modeling  and  Simulations 

We  consider  2-D  breast  models  simulated  in  two  steps.  In 
the  first  step,  the  electromagnetic  field  inside  the  breast  model 
is  simulated  and  the  SAR  distribution  is  calculated  based  on 
the  simulated  electromagnetic  field.  The  second  step  is  for  the 
acoustic  wave  simulation,  where  the  SAR  distribution  obtained 
in  the  first  step  is  used  as  the  acoustic  pressure  source  through 
the  thermal  expansion  coefficient.  In  both  steps,  the  FDTD 
method  is  used  for  the  simulation  examples. 

A.  Electromagnetic  Model  and  Simulation 

For  simulation  purposes,  the  2-D  electromagnetic  breast 
model  used  is  as  shown  in  Fig.  6(a).  The  breast  model  is  a 
10  cm  in  diameter  semicircle,  which  includes  the  skin,  breast 
fatty  tissue,  glandular  tissues,  and  chest  wall  (muscle).  A 
1-mm-diameter  tumor  is  embedded  below  the  skin.  The  di¬ 
electric  properties  of  the  breast  tissues  as  well  as  tumor  at  the 


time  (jis) 

Fig.  7.  Gaussian  modulated  microwave  source. 


microwave  frequency  fi(i  =  1,  •  •  • ,  M)  were  calculated  based 
on  the  Cole-Cole  model  in  (7).  The  dielectric  properties  of  the 
normal  breast  fatty  tissue  are  assumed  random  with  a  variation 
of  ±10%  around  the  nominal  values.  The  dielectric  constants 
of  glandular  tissues  are  between  er  =  11  and  er  =  15. 

Fig.  7  shows  a  Gaussian  modulated  electromagnetic  wave 
used  to  irradiate  the  breast  from  the  top  of  the  model,  as  shown  in 
Fig.  6(a).  The  time  duration  for  the  Gaussian  pulse  is  1  ps.  The 
electromagnetic  field  is  simulated  using  the  FDTD  method  [30], 
[31].  The  grid  cell  size  used  by  FDTD  is  0.5  x  0.5  mm  and  the 
computational  region  is  terminated  by  perfectly  matched  layer 
(PML)  absorbing  boundary  conditions  [44],  [45]. 

The  SAR  distribution  is  given  as  [32],  [33] 


SAR(r) 


^(r)£2(r) 

2p(r) 


(32) 


where  cr(r)  is  the  conductivity  of  the  biological  tissues  at  loca¬ 
tion  r,  E( r)  is  the  electric  field  at  location  r,  and  p(r)  is  the 
mass  density  of  the  biological  tissues  at  location  r. 


B.  Acoustic  Model  and  Simulation 

In  the  microwave-induced  TAI  system,  the  microwave  energy 
is  small,  and  as  a  result,  the  acoustic  pressure  field  induced  by 
microwave  is  also  small.  So,  the  nonlinear  acoustic  effect  does 
not  need  to  be  considered  in  the  TAI  system.  For  example,  it  is 
shown  in  [46]  that  the  temperature  rise  is  about  0.1  mK  and  the 
acoustic  pressure  change  is  only  about  100  Pa  in  the  microwave- 
induced  TAI  system.  The  shock  distance  in  breast  tissues  is  [47] 

B  9  B  2 

v  —  2 A  t  \  _  2 A  m  P0°  t  C 

2(s2  +  1)  P  2(^:  +  1)  P  /max 

(33) 

where  (B/A)(&  10)  is  the  nonlinear  factor  of  the  breast  tis¬ 
sues,  po(~  1000  kg/m3)  is  the  mass  density  of  the  breast  tis¬ 
sues,  and  c(«  1500  m/s)  is  the  sound  speed  inside  the  breast 
tissues  [14].  p  is  the  acoustic  pressure  rise,  and  Amin  and  /max 
are  the  minimal  acoustic  wavelength  and  the  maximal  acoustic 
frequency  of  the  thermal  acoustic  signal,  respectively.  For  our 
breast  model,  the  acoustic  pressure  rise  is  p  —  100  Pa,  and  the 
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TABLE  II 

Acoustic  Parameters  for  Biological  Tissues. 

(*  /  Is  the  Acoustic  Frequency  and  the  Unit  Is  Megahertz) 


Tissue 

Breast 

Skin 

Muscle 

Tumor 

p  {kg/m6) 

1020 

1100 

1041 

1041 

c  {m/s) 

1510 

1537 

1580 

1580 

a *  {dB/cm) 

0.75/lb 

3.5 

0.57/ 

0.57/ 

P  (1/°C) 

3E-4 

3E-4 

3E-4 

3E-4 

Cp  ( J/(°C  ■  kg)) 

3550 

3500 

3510 

3510 

maximal  acoustic  frequency  is  /max  =  500  KHz.  By  substi¬ 
tuting  the  parameters  into  (33),  we  obtain  the  shock  distance  in 
breast  tissues  to  be 

_  2A  m  P0°2  m  C 
2(^  +  l)  V  /max 

5  1000  •  15002  1500 

”  2(5  +  1)  100  500  x  103 

=  2.8  x  104  m.  (34) 

Because  the  size  of  our  breast  model  is  only  10  cm,  which  is 
much  smaller  than  the  shock  distance,  it  is  reasonable  to  ig¬ 
nore  the  nonlinear  acoustic  effect  in  the  microwave-induced  TAI 
system. 

The  two  basic  linear  acoustic  wave  generation  equations  are 
[13] 

P^u(r,  t)  =  —  Vp(r,  t)  (35) 

and 

\  d  d 

V  •  u(r, i)  =  --^—p(r,t)  +  ap(r,t)+/3—T(r,t)  (36) 

where  u(r,t)  is  the  acoustic  velocity  vector,  p(r,t)  is  the 
acoustic  pressure  field,  p  is  the  mass  density,  a  is  the  at¬ 
tenuation  coefficient,  6  is  the  thermal  expansion  coefficient, 
and  T(r,  t)  is  the  temperature.  The  values  for  these  acoustic 
properties  for  different  breast  tissues  are  listed  in  Table  II  [14], 
[46],  [48]-[50].  The  attenuation  coefficient  is  calculated  with 
/  =  0.15  MHz.  The  values  for  the  tumor  are  approximated 
using  the  values  for  muscle  because  we  cannot  find  the  values 
specific  to  the  tumor. 

Because  the  duration  of  the  microwave  pulse  is  much  shorter 
than  the  thermal  diffusion  time,  thermal  diffusion  can  be  ne¬ 
glected  [13],  and  the  thermal  equation  is 

CpJ^T(r,  i)  =  SAR(r,  t)  (37) 

where  Cp  is  the  specific  heat.  Substituting  (37)  into  (36)  gives 

\  d  3 

V-u(r,£)  = - 2— p(r,t)  +  ap(r,t)  +  -^-SAR(r,t ).  (38) 

pc  ot  Up 

FDTD  is  used  again  to  compute  the  thermal  acoustic  wave 
based  on  (35)  and  (38).  More  details  about  FDTD  for  acoustic 
simulations  can  be  found  in  [34],  [35],  and  [5 1]— [57]. 

The  breast  model  for  the  acoustic  simulation  is  constructed 
similarly  to  the  model  for  electromagnetic  simulation.  The 


Thermal  Acoustic  Signals  from  Tumor 


Normalized  Spectrums  of  The  Thermal  Acoustic  Signals 


Frequency  (KHz) 

(b) 

Fig.  8.  Thermal  acoustic  signals  at  different  stimulating  frequencies  /  =  200, 
400,  600,  and  800  MHz.  (a)  Thermal  acoustic  responses  from  tumor  only  and 
(b)  the  normalized  spectrums  of  the  signals  in  (a). 


velocities  of  the  normal  fatty  breast  tissue  are  also  assumed 
random  with  a  variation  of  ±5%  around  average  values,  as 
shown  in  Fig.  6(b).  An  acoustic  sensor  array  with  35  elements 
deployed  uniformly  around  the  breast  model  is  used  to  record 
the  thermal  acoustic  signals.  The  distance  between  neighboring 
acoustic  sensors  is  4  mm.  The  grid  cell  size  used  by  the  acoustic 
FDTD  is  0.1  x  0.1  mm  and  the  computational  region  is  ter¬ 
minated  by  PML  absorbing  boundary  conditions  [55]-[57]. 
Note  that  the  size  of  the  FDTD  cell  for  acoustic  simulation 
is  much  finer  than  that  of  the  FDTD  cell  for  electromagnetic 
simulation  because  the  wavelength  of  an  acoustic  wave  is  much 
smaller  than  that  of  a  microwave.  The  SAR  distribution  data 
is  interpolated  to  achieve  the  designed  grid  resolution  for  the 
acoustic  breast  model. 

V.  Numerical  Examples 

At  the  beginning  of  this  section,  the  typical  microwave-in¬ 
duced  thermal  acoustic  responses  from  the  tumor  are  plotted 
in  Fig.  8(a)  for  stimulating  frequencies  of  /  =  200,  400,  600, 
and  800  MHz.  The  signals  are  simulated  based  on  the  afore¬ 
mentioned  2-D  breast  model.  To  obtain  the  signals,  we  perform 
the  simulation  twice  at  each  stimulating  frequency,  with  and 
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Fig.  9.  Imaging  results  for  the  case  of  a  single  1  mm-diameter  tumor, 
(a)  MART,  (b)  DAS,  (c)  SART  at  stimulating  frequency  /  =  200  MHz, 
(d)  SART  at  stimulating  frequency  /  =  400  MHz,  (e)  SART  at  stimulating  fre¬ 
quency  /  =  600  MHz,  and  (f)  SART  at  stimulating  frequency  /  =  800  MHz. 


Fig.  10.  Imaging  results  for  the  two  1.5 -mm-diameter  tumors  case,  (a)  Breast 
model,  (b)  MART,  (c)  DAS,  (d)  zoom  in  of  (b),  (e)  SART  at  stimulating  fre¬ 
quency  /  =  300  MHz,  and  (f)  SART  at  stimulating  frequency  /  =  700  MHz. 


without  the  tumor,  and  record  the  acoustic  signals  in  an  acoustic 
sensor.  The  difference  of  the  two  received  signals  is  referred  to 
as  the  thermal  acoustic  response  only  from  the  tumor  at  the  stim¬ 
ulating  frequency.  It  can  be  seen  that  the  thermal  acoustic  re¬ 
sponses  from  the  tumor  at  different  stimulating  frequencies  are 
similar  to  one  another.  The  figure  also  shows  that  the  thermal 
acoustic  signals  are  wideband  bipolar  pulses,  with  a  large  pos¬ 
itive  peak  and  a  large  negative  peak.  Fig.  8(b)  shows  the  nor¬ 
malized  spectra  of  the  acoustic  signals  corresponding  to  the  ex¬ 
citations  in  Fig.  8(a).  It  is  seen  that  the  frequency  range  of  the 
acoustic  signals  is  about  from  1  to  400  KHz.  The  dominant  band 
(3-dB  band)  of  the  signals  ranges  from  10  to  180  KHz,  and  the 
corresponding  acoustic  wavelength  ranges  from  150  to  8  mm  in 
the  breast  tissues. 

Several  numerical  examples  are  used  in  this  section  to  demon¬ 
strate  the  effectiveness  of  MART.  The  thermal  acoustic  sig¬ 
nals  are  simulated  based  on  the  2-D  breast  model  with  tumor 
only.  For  comparison  purposes,  the  DAS  method  is  applied  to 
the  same  data  set  also.  We  also  present  the  imaging  results  for 
the  single-frequency  microwave-induced  TAI  at  different  stim¬ 
ulating  frequencies.  The  corresponding  image  reconstruction 
method  is  referred  to  as  the  single-frequency  adaptive  and  robust 
technique  (SART),  which  is  similar  with  MART  but  without 
Stage  II  of  MART.  In  SART,  the  RCB  is  used  to  estimate  the 
thermal  acoustic  waveform  at  a  certain  stimulating  frequency 
just  like  in  Stage  I  of  MART.  Then  the  peak  search  method  used 
in  MART  Stage  III  is  applied  to  the  estimated  waveform  to  de¬ 
termine  the  image  intensity. 


In  the  first  example,  a  1 -mm-diameter  tumor  is  embedded  in 
the  breast  model  at  the  location  (x  =  70  mm,  y  =  60  mm). 
This  is  the  challenging  case  of  early  breast  cancer  detection 
because  of  the  small  tumor  size.  Fig.  9(a)  and  (b)  shows  the 
imaging  results  for  MART  and  DAS,  respectively.  The  tumor 
is  shown  clearly  in  the  MART  image  [Fig.  9(a)],  and  the  size 
and  location  of  the  tumor  is  accurate.  Because  of  the  high  side- 
lobe,  poor  resolution,  and  poor  interference  rejection  capability 
of  the  DAS  method,  the  tumor  is  essentially  missed  by  DAS  as 
shown  in  Fig.  9(b).  Fig.  9(c)-(f)  shows  the  imaging  results  for 
SART  at  the  stimulating  frequencies  /  =  200,  400,  600,  and 
800  MHz,  respectively.  The  figures  show  that  SART  can  de¬ 
termine  the  tumor  correctly,  but  some  clutter  shows  up  in  the 
SART  images.  Note  that  the  clutter  shows  up  at  different  loca¬ 
tions  with  different  stimulating  frequencies.  By  comparing  the 
images  for  MART  and  SART,  it  can  be  seen  that  the  clutter  are 
effectively  suppressed  by  MART  when  multiple  stimulating  fre¬ 
quencies  are  used. 

In  the  second  numerical  example,  two  small  1.5-mm-di- 
ameter  tumors  are  set  inside  the  breast  model  as  shown  in 
Fig.  10(a).  Their  locations  are  at  (x  =  70  mm,  y  =  60  mm)  and 
(x  =  75  mm,  y  =  62.5  mm).  The  distance  between  the  two 
tumors  is  4  mm.  The  imaging  results  using  MART  and  DAS 
are  shown  in  Fig.  10(b)  and  (c),  respectively.  The  two  tumors 
are  seen  clearly  in  the  MART  image.  To  show  them  clearly 
we  zoom  in  onto  the  tumor  locations  in  Fig.  10(d),  where  the 
two  black  circles  mark  the  actual  sizes  and  locations  of  the 
two  tumors.  It  is  shown  that  MART  can  be  used  to  determine 
the  locations  and  sizes  of  the  two  tumors  accurately.  The  DAS 
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image  contains  much  clutter.  The  two  tumors  cannot  be  sepa¬ 
rated  clearly  in  the  DAS  image  because  of  the  poor  resolution 
of  DAS.  Fig.  10(e)  and  (f)  shows  the  imaging  results  of  SART 
at  stimulation  frequencies  /  =  300  and  700  MHz,  respectively. 
The  tumors  can  be  seen  in  both  of  the  SART  images,  but  clutter 
shows  up  between  the  two  tumors  in  Fig.  10(e)  and  (f),  and  the 
sizes  of  the  two  tumors  in  Fig.  10(f)  are  larger  than  their  true 
sizes. 

VI.  Conclusion 

An  investigation  of  using  a  multifrequency  microwave-in¬ 
duced  TAI  system  for  early  breast  cancer  detection  has  been 
reported  in  this  paper.  The  frequency  band  for  this  system  has 
been  given  based  on  the  cutoff  frequency  of  the  human  breast. 
A  simplified  semicircular  dielectric  waveguide  mode  was  used 
to  calculate  the  cutoff  frequency  in  this  paper.  By  studying  the 
microwave  energy  absorption  properties  of  breast  tissue  and 
tumor,  we  have  shown  that  the  multifrequency  microwave-in¬ 
duced  TAI  can  offer  higher  SNR,  higher  imaging  contrast,  and 
more  effective  clutter  suppression  capability  than  the  traditional 
single-frequency  microwave-induced  TAI.  A  MART  has  been 
presented  for  image  formation.  This  data-adaptive  algorithm 
can  achieve  better  resolution  and  better  interference  rejection 
capability  than  its  data-independent  counterparts,  such  as  DAS. 
The  feasibility  of  this  multifrequency  microwave-induced  TAI 
system  as  well  as  the  performance  of  the  proposed  image  re¬ 
construction  algorithm  for  early  breast  cancer  detection  have 
been  demonstrated  by  using  2-D  numerical  electromagnetic  and 
acoustic  breast  models.  The  absorbed  microwave  energy  and  the 
thermal  acoustic  field  in  the  breast  models  have  been  simulated 
using  the  FDTD  method.  Numerical  examples  have  been  used 
to  demonstrate  the  excellent  performance  of  this  multifrequency 
technique.  More  advanced  models  are  being  developed  to  fur¬ 
ther  investigate  and  validate  the  preferential  imaging  capability 
of  the  technique  for  early  breast  cancer  detection. 

References 

[1]  E.  C.  Fear,  P.  M.  Meaney,  and  M.  A.  Stuchly,  “Microwaves  for  breast 
cancer  detection,”  IEEE  Potentials,  vol.  22,  no.  1,  pp.  12-18,  Feb./Mar. 
2003. 

[2]  P.  M.  Meaney,  M.  W.  Fanning,  D.  Li,  S.  P.  Poplack,  and  K.  D.  Paulsen, 
“A  clinical  prototype  for  active  microwave  imaging  of  the  breast,” 
IEEE  Trans.  Microw.  Theory  Tech.,  vol.  48,  no.  11,  pp.  1841-1853, 
Nov.  2000. 

[3]  E.  J.  Bond,  X.  Li,  S.  C.  Hagness,  and  B.  D.  Van  Veen,  “Microwave 
imaging  via  space-time  beamforming  for  early  detection  of  breast 
cancer,”  IEEE  Trans.  Antennas  Propagat.,  vol.  51,  no.  8,  pp. 
1690-1705,  Aug.  2003. 

[4]  B.  Guo,  Y.  Wang,  J.  Li,  P.  Stoica,  and  R.  Wu,  “Microwave  imaging  via 
adaptive  beamforming  methods  for  breast  cancer  detection,”  J.  Elec- 
tromagn.  Waves  Appl.,  vol.  20,  no.  1,  pp.  53-63,  Jan.  2006. 

[5]  Y.  Xie,  B.  Guo,  L.  Xu,  J.  Li,  and  P.  Stoica,  “Multi-static  adaptive 
microwave  imaging  for  early  breast  cancer  detection,”  IEEE  Trans. 
Biomed.  Eng.,  vol.  53,  no.  8,  pp.  1647-1657,  Aug.  2006. 

[6]  D.  Li,  P.  M.  Meaney,  and  K.  D.  Paulsen,  “Conformal  microwave 
imaging  for  breast  cancer  detection,”  IEEE  Trans.  Microw.  Theory 
Tech.,  vol.  51,  no.  4,  pp.  1179-1186,  Apr.  2003. 

[7]  Q.  F.  Fang,  P.  M.  Meaney,  S.  D.  Geimer,  A.  V.  Streltsov,  and  K.  D. 
Paulsen,  “Microwave  image  reconstruction  from  3-D  fields  coupled  to 
2-D  parameter  estimation,”  IEEE  Trans.  Med.  Imag.,  vol.  23,  no.  4,  pp. 
475-484,  Apr.  2004. 

[8]  Q.  H.  Liu,  Z.  Q.  Zhang,  T.  T.  Wang,  J.  A.  Bryan,  G.  A.  Ybarra,  L.  W. 
Nolte,  and  W.  T.  Joines,  “Active  microwave  imaging  I — 2-D  forward 
and  inverse  scattering  methods,”  IEEE  Trans.  Microw.  Theory  Tech., 
vol.  50,  no.  1,  pp.  123-133,  Jan.  2002. 


[9]  Z.  Q.  Zhang,  Q.  H.  Liu,  C.  Xiao,  E.  Ward,  G.  Ybarra,  and  W.  T.  Joines, 
“Microwave  breast  imaging:  3-D  forward  scattering  simulation,”  IEEE 
Trans.  Biomed.  Eng.,  vol.  50,  no.  10,  pp.  1180-1189,  Oct.  2003. 

[10]  R.  K.  Mishra,  A.  Swarup,  and  J.  M.  Thomas,  “Dielectric  properties 
of  normal  and  malignant  human  breast  tissues  at  radiowave  and  mi¬ 
crowave  frequencies,”  Ind.  J.  Biochem.  Biophys.,  vol.  21,  pp.  76-79, 
Feb.  1984. 

[11]  A.  J.  Surowiec,  S.  S.  Stuchly,  J.  R.  Barr,  and  A.  Swarup,  “Dielec¬ 
tric  properties  of  breast  carcinoma  and  the  surrounding  tissues,”  IEEE 
Trans.  Biomed.  Eng.,  vol.  35,  no.  4,  pp.  257-263,  Apr.  1988. 

[12]  W.  T.  Joines,  Y.  Zhang,  C.  X.  Li,  and  R.  L.  Jirtel,  “The  measured  elec¬ 
trical  properties  of  normal  and  malignant  human  tissues  from  50  to  900 
MHz,”  Med.  Phys.,  vol.  21,  no.  4,  pp.  547-550,  Apr.  1994. 

[13]  M.  Xu  and  L.  V.  Wang,  “Time-domain  reconstruction  for  thermoa¬ 
coustic  tomography  in  a  speherical  geometry,”  IEEE  Trans.  Med. 
Imag.,  vol.  21,  no.  7,  pp.  814-822,  Jul.  2002. 

[14]  T.  Szabo,  Diagnostic  Ultrasound  Imaging:  Inside  Out.  Burlington, 
MA:  Elsevier,  2004. 

[15]  F.  W.  Kremkau,  Diagnostic  Ultrasound:  Principles  and  Instruments . 
Philadelphia,  PA:  Saunders,  1993. 

[16]  E.  Steen  and  B.  Olstad,  “Volume  rendering  of  3-D  medical  ultrasound 
data  using  direct  feature  mapping,”  IEEE  Trans.  Med.  Imag.,  vol.  13, 
no.  6,  pp.  517-525,  Jun.  1994. 

[17]  W.  Joines,  R.  Jirtle,  M.  Rafal,  and  D.  Schaeffer,  “Microwave  power  ab¬ 
sorption  diffrences  between  normal  and  malignant  tissue,”  Int.  J.  Ra- 
diat.  Oncol.  Biol.  Phys.,  vol.  6,  no.  6,  pp.  681-687,  1980. 

[18]  Y.  Xu  and  L.  H.  Wang,  “Effects  of  acoustic  heterogeneity  in  breast 
thermoacoustic  tomography,”  IEEE  Trans.  Ultrason.  Ferroelect.  Freq. 
Control,  vol.  50,  no.  9,  pp.  1134-1146,  Sep.  2003. 

[19]  M.  Xu,  Y.  Xu,  and  L.  V.  Wang,  “Time-domain  reconstruction  algo¬ 
rithms  and  numerical  simulations  for  thermoacoustic  tomography  in 
various  geometries,”  IEEE  Trans.  Biomed.  Eng.,  vol.  50,  no.  9,  pp. 
1086-1099,  Sep.  2003. 

[20]  Y.  Xu,  P.  Kuchment,  and  L.  H.  Wang,  “Limited-view  thermoacoustic 
tomography  and  reconstruction  by  truncated  conjugate  gradient,”  Med. 
Phys.,  vol.  31,  pp.  724-733,  Apr.  2004. 

[21]  R.  A.  Kruger,  K.  K.  Kopecky,  A.  M.  Aisen,  D.  R.  Reinecke,  G.  A. 
Kruger,  and  W.  L.  Kiser,  “Thermoacoustic  CT  with  radio  waves:  A 
medical  imaging  paradigm,”  Radiology,  vol.  211,  pp.  275-278,  Apr. 
1999. 

[22]  R.  A.  Kruger,  K.  D.  Miller,  H.  E.  Reynolds,  W.  L.  Kiser,  D.  R.  Rei¬ 
necke,  and  G.  A.  Kruger,  “Breast  cancer  in  vivo:  Contrast  enhancement 
with  thermoacoustic  CT  at  434  MHz — Feasibility  study,”  Radiology, 
vol.  216,  pp.  279-283,  Jul.  2000. 

[23]  R.  A.  Kruger,  W.  L.  Kiser,  K.  D.  Miller,  H.  E.  Reynolds,  D.  R.  Rei¬ 
necke,  and  G.  A.  Kruger,  “Thermoacoustic  CT,”  in  Proc.  IEEE  MTT-S 
Int.  Microwave  Symp.  Dig.,  Jun.  2000,  vol.  2,  pp.  11-16. 

[24]  Y.  Xu,  D.  Feng,  and  L.  H.  Wang,  “Exact  frequency-domain  reconstruc¬ 
tion  for  thermoacoustic  tomography  I:  Planar  geometry,”  IEEE  Trans. 
Med.  Imag.,  vol.  21,  no.  7,  pp.  823-828,  Jul.  2002. 

[25]  Y.  Xu,  M.  Xu,  and  L.  H.  Wang,  “Exact  frequency-domain  reconstruc¬ 
tion  for  thermoacoustic  tomography  II:  Cylindrical  geometry,”  IEEE 
Trans.  Med.  Imag.,  vol.  21,  no.  7,  pp.  829-833,  Jul.  2002. 

[26]  J.  Li,  P.  Stoica,  and  Z.  Wang,  “On  robust  capon  beamforming  and 
diagonal  loading,”  IEEE  Trans.  Signal  Process.,  vol.  51,  no.  7,  pp. 
1702-1715,  Jul.  2003. 

[27]  J.  Li  and  P.  Stoica,  Eds.,  Robust  Adaptive  Beamforming.  New  York: 
Wiley,  2005. 

[28]  Y.  Xie,  B.  Guo,  J.  Li,  and  P.  Stoica,  “Novel  multistatic  adaptive  mi¬ 
crowave  imaging  methods  for  early  breast  cancer  detection,”  EURASIP 
J.  Appl.  Signal  Process.,  vol.  2006,  pp.  1-12,  2006. 

[29]  Y.  Xie,  B.  Guo,  J.  Li,  G.  Ku,  and  L.  V.  Wang,  “Adaptive  and  robutst 
techniques  (ART)  for  thermoacoustic  tomography,”  presented  at  the 
Proc.  40th  Asilomar  Conf.  Signals,  Systems  and  Computers,  Pacific 
Grove,  CA,  Nov.  2006. 

[30]  A.  Taflove  and  S.  C.  Hagness,  Computational  Electrodynamics:  The  Fi¬ 
nite-Difference  Time-Domain  Method,  2nd  ed.  Boston,  MA:  Artech 
House,  2000. 

[31]  D.  M.  Sullivan,  Electromagnetic  Simulation  Using  the  FDTD 
Method.  New  York:  IEEE,  2000. 

[32]  P.  Bernardi,  M.  Cavagnaro,  S.  Pisa,  and  E.  Piuzzi,  “Specific  ab¬ 
sorption  rate  and  temperature  elevation  in  a  subject  exposed  in  the 
far-field  of  radio-frequency  sources  operating  in  the  10-900-MHz 
range,”  IEEE  Trans.  Biomed.  Eng.,  vol.  50,  no.  3,  pp.  295-304,  Mar. 
2003. 

[33]  P.  Bernardi,  M.  Cavagnaro,  S.  Pisa,  and  E.  Piuzzi,  “SAR  distribution 
and  temperature  increase  in  an  anatomical  model  of  the  human  eye 
exposed  to  the  field  radiated  by  the  user  antenna  in  a  wireless  LAN,” 
IEEE  Trans.  Microw.  Theory  Tech.,  vol.  46,  no.  12,  pp.  2074-2082, 
Dec.  1998. 


2010 


IEEE  TRANSACTIONS  ON  BIOMEDICAL  ENGINEERING,  VOL.  54,  NO.  11,  NOVEMBER  2007 


[34]  J.  G.  Maloney  and  K.  E.  Cummings,  “Adaptation  of  FDTD  techniques 
to  acoustic  modeling,”  in  Proc.  11th  Annu.  Rev.  Progress  in  Applied 
Computational  Electromagnetics,  Monterey,  CA,  1995,  vol.  2,  pp. 
724-731. 

[35]  X.  Yuan,  D.  Borup,  J.  Wiskin,  M.  Berggren,  and  S.  A.  Johnson,  “Simu¬ 
lation  of  acoustic  wave  propagation  in  dispersive  media  with  relaxation 
losses  by  using  FDTD  method  with  PML  absorbing  boundary  condi¬ 
tion,”  IEEE  Trans.  Ultrason.,  Ferroelect.,  Freq.  Control,  vol.  46,  no.  1, 
pp.  14-23,  Jan.  1999. 

[36]  C.  A.  Balanis,  Advanced  Engineering  Electromagnetics.  New  York: 
Wiley,  1989. 

[37]  S.  Yang,  H.  Kim,  and  H.  J.  Lee,  “Circular-harmonic  vector  analysis  of 
dielectric  waveguide  with  a  cross-cut-circle  cross  section,”  Appl.  Opt., 
vol.  34,  pp.  7705-7713,  Nov.  1995. 

[38]  S.  Gabriel,  R.  W.  Lau,  and  C.  Gabriel,  “The  dielectric  properties  of 
biological  tissues:  III.  Parametric  models  for  the  dielectric  spectrum  of 
tissues,”  Phys.  Med.  Biol.,  vol.  41,  pp.  2271-2293,  Nov.  1996. 

[39]  C.  Gabriel  and  S.  Gabriel,  “Compilation  of  the  dielectric  proper¬ 
ties  of  body  tissues  at  rf  and  microwave  frequencies,”  [Online]. 
Available:  http://www.brooks.af.mil/AFRL/HED/hedr/reports/dielec- 
tric/home.html 

[40]  “Dielectric  properties  of  body  tissues,”  [Online].  Available: 
http://niremf.ifac.cnr.it/tissprop/ 

[41]  X.  Li  and  S.  C.  Hagness,  “A  confocal  microwave  imaging  algorithm  for 
breast  cancer  detection,”  IEEE  Microw.  Wireless  Compon.  Lett.,  vol. 
11,  no.  3,  pp.  130-132,  Mar.  2001. 

[42]  P.  Stoica  and  R.  L.  Moses,  Spectral  Analysis  of  Signals.  Upper  Saddle 
River,  NJ:  Prentice-Hall,  2005. 

[43]  C.  G.  A.  Hoelen  and  F.  F.  M.  de  Mul,  “Image  reconstruction  for 
photoacoustic  scanning  of  tissue  structures,”  Appl.  Opt.,  vol.  39,  pp. 
5872-5883,  Nov.  2000. 

[44]  S.  D.  Gedney,  “An  anisotropic  perfectly  matched  layer-absorbing 
medium  for  the  truncation  of  FDTD  lattices,”  IEEE  Trans.  Antennas 
Propagat.,  vol.  44,  no.  12,  pp.  1630-1639,  Dec.  1996. 

[45]  Z.  S.  Sacks,  D.  M.  Kingsland,  R.  Lee,  and  J.  F.  Lee,  “A  perfectly 
matched  anisotropic  absorber  for  use  as  an  absorbing  boundary  condi¬ 
tion,”  IEEE  Trans.  Antennas  Propagat.,  vol.  43,  no.  12,  pp.  1460-1463, 
Dec.  1995. 

[46]  M.  Xu  and  L.  H.  Wang,  “Pulsed-microwave-induced  thermoacoustic 
tomography:  Filtered  backprojection  in  a  circular  measurement  con¬ 
figuration,”  Med.  Phys.,  vol.  29,  no.  8,  pp.  1661-1669,  Aug.  2002. 

[47]  M.  F.  Hamilton  and  D.  T.  Blackstock,  Nonlinear  Acoustics.  San 
Diego,  CA:  Acadmic,  1997. 

[48]  F.  A.  Duck,  Physical  Properties  of  Tissue.  London,  U.K.:  Academic, 
1990. 

[49]  M.  L.  Palmeri  and  K.  R.  Nightingale,  “On  the  thermal  effects  associated 
with  radiation  force  imaging  of  soft  tissue,”  IEEE  Trans.  Ultrason., 
Ferroelect.,  Freq.  Control,  vol.  51,  no.  5,  pp.  551-565,  May  2004. 

[50]  M.  Converse,  E.  J.  Bond,  S.  C.  Hagness,  and  B.  D.  Van  Veen,  “Ultra- 
wide-band  microwave  space-time  beamforming  for  hyperthermia  treat¬ 
ment  of  breast  cancer:  A  computational  feasibility  study,”  IEEE  Trans. 
Microw.  Theory  Tech.,  vol.  52,  no.  8,  pp.  1876-1889,  Aug.  2004. 

[51]  C.  K.  W.  Tam  and  J.  C.  Webb,  “Dispersion-relation-preserving  finite- 
difference  schemes  for  computational  acoustics,”  J.  Comput.  Phys. ,  vol. 
107,  no.  8,  pp.  262-281,  Aug.  1993. 

[52]  S.  V.  Tsynkov,  “Numerical  solution  of  problems  on  unbounded  do¬ 
mains.  A  review,”  Appl.  Numer.  Math.,  vol.  27,  pp.  465-532,  Aug. 
1998. 

[53]  T.  Colonius,  S.  K.  Lele,  and  P.  Moin,  “Boundary- conditions  for  direct 
computation  of  aerodynamic  sound  generation,”  AIAA  J.,  vol.  31,  pp. 
1574-1582,  Sep.  1993. 

[54]  C.  W.  Rowley  and  T.  Colonius,  “Discretely  nonreflecting  boundary 
conditions  for  linear  hyperbolic  systems,”  J.  Comput.  Phys.,  vol.  157, 
pp.  500-538,  Jan.  2000. 

[55]  Q.  H.  Liu  and  J.  Tao,  “The  perfectly  matched  layer  for  acoustic  waves 
in  absorptive  media,”  J.  Acoust.  Soc.  Amer.,  vol.  102,  pp.  2072-2082, 
October  1997. 

[56]  X.  Yuan,  D.  Borup,  J.  W.  Wiskin,  M.  Berggren,  R.  Eidens,  and  S.  A. 
Johnson,  “Formulation  and  validation  of  Berenger’s  PML  absorbing 
boundary  for  the  FDTD  simulation  of  acoustic  scattering,”  IEEE  Trans. 
Ultrason.,  Ferroelect.,  Freq.  Control,  vol.  44,  no.  7,  pp.  816-822,  Jul. 
1997. 

[57]  T.  K.  Katsibas  and  C.  S.  Antonopoulos,  “A  general  form  of  perfectly 
matched  layers  for  three-dimensional  problems  of  acoustic  scattering 
in  lossless  and  lossy  fluid  media,”  IEEE  Trans.  Ultrason.,  Ferroelect., 
Freq.  Control,  vol.  51,  no.  8,  pp.  964-972,  Aug.  2004. 


Bin  Guo  (S’06-M’07)  received  the  B.E.  and  M.Sc. 
degrees  in  electrical  engineering  from  Xian  Jiaotong 
University,  Xian,  China,  in  1997  and  2000,  respec¬ 
tively.  He  is  currently  pursuing  the  Ph.D.  degree  in 
electrical  engineering  in  the  Department  of  Electrical 
and  Computer  Engineering,  University  of  Florida, 
Gainesville. 

From  April  2002  to  July  2003,  he  was  an  Asso¬ 
ciate  Research  Scientist  with  the  Temasek  Labora¬ 
tories,  National  University  of  Singapore,  Singapore. 
Since  August  2003,  he  has  been  a  Research  Assis¬ 
tant  with  the  Department  of  Electrical  and  Computer  Engineering,  University  of 
Florida.  His  current  research  interests  include  biomedical  applications  of  signal 
processing,  microwave  imaging,  and  computational  electromagnetics. 


Jian  Li  (S’SS-M^O-SM^-F’CU)  received  the 
M.Sc.  and  Ph.D.  degrees  in  electrical  engineering 
from  The  Ohio  State  University,  Columbus,  in  1987 
and  1991,  respectively. 

From  July  1991  to  June  1993,  she  was  an  Assis¬ 
tant  Professor  with  the  Department  of  Electrical  En¬ 
gineering,  University  of  Kentucky,  Lexington.  Since 
August  1993,  she  has  been  with  the  Department  of 
Electrical  and  Computer  Engineering,  University  of 
Florida,  Gainesville,  where  she  is  currently  a  Pro¬ 
fessor.  Her  current  research  interests  include  spectral 
estimation,  statistical  and  array  signal  processing,  and  their  applications. 

Dr.  Li  is  a  fellow  of  Institution  of  Electrical  Engineers  (IEE).  She  received 
the  1994  National  Science  Foundation  Young  Investigator  Award  and  the  1996 
Office  of  Naval  Research  Young  Investigator  Award.  She  has  been  a  member  of 
the  Editorial  Board  of  Signal  Processing,  a  publication  of  the  European  Associ¬ 
ation  for  Signal  Processing  (EURASIP)  since  2005.  She  is  presently  a  member 
of  two  of  the  IEEE  Signal  Processing  Society  technical  committees:  the  Signal 
Processing  Theory  and  Methods  (SPTM)  Technical  Committee  and  the  Sensor 
Array  and  Multichannel  (SAM)  Technical  Committee. 


Henry  Zmuda  (M’75)  received  the  B.E.  degree 
from  the  Stevens  Institute  of  Technology,  Hoboken, 
NJ,  in  1979,  and  the  M.S.  and  Ph.D.  degrees  from 
Cornell  University,  Ithaca,  NY,  in  1982  and  1984, 
respectively,  all  in  electrical  engineering. 

From  1984  to  1995,  he  was  a  member  of  the 
faculty  at  Stevens  Institute  of  Technology.  In  1995, 
he  joined  the  Electrical  Engineering  Faculty,  Univer¬ 
sity  of  Florida.  His  research  interests  are  primarily 
in  the  area  of  microwave  photonics  and  applied 
electromagnetics. 


Mark  Sheplak  received  the  B.S.  degree  in  1989,  the 
M.S.  degree  in  1992,  and  the  Ph.D.  degree  in  1995 
in  mechanical  engineering  from  Syracuse  University, 
Syracuse,  NY. 

He  is  an  Associate  Professor  in  the  Department  of 
Aerospace  and  Mechanical  Engineering  and  an  Af¬ 
filiate  Associate  Professor  of  electrical  and  computer 
engineering  at  the  University  of  Florida  (UF).  Prior 
to  joining  UF  in  1998,  he  was  a  Postdoctoral  Asso¬ 
ciate  at  the  Massachusetts  Institute  of  Technology’s 
Microsystems  Technology  Laboratories,  Cambridge, 
from  1995  to  1998.  During  his  Ph.D.  studies,  he  was  a  GSRP  Fellow  at  the 
NASA-Langley  Research  Center,  Hampton,  VA,  from  1992  to  1995.  His  current 
research  focuses  on  the  design,  fabrication,  and  characterization  of  high-perfor¬ 
mance,  instrumentation-grade,  MEMS-based  sensors  and  actuators  that  enable 
the  measurement,  modeling,  and  control  of  various  physical  properties. 


Adaptive  and  Robust  Techniques  (ART)  for 
Thermoacoustic  Tomography  in  Breast  Cancer 

Detection 


Yao  Xie*  Bin  Guo*  Jian  Li*  Geng  Ku*  Lihong  Wang* 

*Dept.  of  Electrical  and  Computer  Engineering  *Dept.  of  Biomedical  Engineering 

University  of  Florida  Texas  A&M  University 

P.  O.  Box  116130,  Gainesville,  FL  32611,  USA.  3120  TAMU,  College  Station,  TX  77843-3120,  USA. 


1  Abstract — We  present  Adaptive  and  Robust  Techniques 
(ART)  for  thermoacoustic  tomography  with  application  in  breast 
cancer  detection.  TAT  combines  high  contrast  and  high  resolu¬ 
tion.  The  current  image  reconstruction  methods  used  for  TAT, 
such  as  the  widely  used  Delay-and-Sum  (DAS)  approach,  are 
data-independent  and  suffer  from  low  resolution,  high  sidelobe 
levels,  and  poor  interference  rejection  capabilities.  The  data- 
adaptive  ART  have  much  better  performance  than  their  data- 
independent  counterparts.  By  allowing  certain  uncertainties, 
ART  can  mitigate  the  amplitude  and  phase  distortion  problems 
encountered  in  TAT.  The  excellent  performance  of  ART  is 
demonstrated  using  both  simulated  and  experimentally  measured 
data. 

I.  Introduction 

Thermoacoustic  tomography  (TAT)  [1]  is  a  new  technology 
with  great  promise  in  a  wide  span  of  biomedical  applications. 
In  TAT,  an  image  of  the  tissue  radiation  absorption  properties 
is  reconstructed  from  the  recorded  thermoacoustic  signals. 
Such  an  image  reveals  the  physiological  and  pathological 
status  of  the  tissue,  which  has  been  used  in  many  applications 
including  breast  cancer  detection.  Microwave-induced  TAT 
combines  the  high  contrast  of  biological  tissue  in  electro¬ 
magnetic  frequency  band  and  the  millimeter  range  spatial 
resolution  due  to  the  acoustic  signals  used  [2]. 

Developing  accurate  and  robust  image  reconstruction  meth¬ 
ods  is  one  of  the  key  challenges  encountered  in  TAT.  Time- 
domain  approximate  reconstruction  algorithms  such  as  the 
DAS  [3]  (weighted  or  unweighted)  type  of  data-independent 
approaches  are  widely  used.  They  need  little  prior  information 
on  the  tissue  for  image  reconstruction  and  can  be  fast  and 
simple  to  implement  to  process  the  wideband  acoustic  sig¬ 
nals.  Although  not  based  on  the  exact  solution,  they  provide 
similar  image  qualities  to  those  of  the  exact  reconstruction 
algorithms.  However,  these  data-independent  methods  suffer 
from  poor  resolution  and  high  sidelobe  level  problems.  Data- 
adaptive  approaches,  such  as  the  recently  introduced  Robust 
Capon  Beamforming  (RCB)  method  [4],  can  have  much  better 
resolution  and  much  better  interference  rejection  capability 
than  their  data-independent  counterparts. 
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A  common  assumption  of  the  existing  methods  is  that  the 
surrounding  tissue  is  acoustically  homogeneous.  This  approx¬ 
imation  is  inadequate  in  many  medical  imaging  applications. 
According  to  previous  studies,  the  sound  speed  in  human 
female  breast  varies  widely  from  1430  m/s  to  1570  m/s  around 
the  commonly  assumed  speed  of  1510  m/s  [5].  The  heteroge¬ 
neous  acoustic  properties  of  biological  tissues  cause  amplitude 
and  phase  distortions  in  the  recorded  acoustic  signals,  which 
can  result  in  significant  degradations  in  image  quality  [6].  Al¬ 
though  various  methods  have  been  proposed  for  the  amplitude 
and  phase  distortions  in  Ultra-sound  Tomography  (UT),  the 
problem  in  TAT  is  somewhat  different:  in  TAT,  for  the  breast 
tissue  with  a  relatively  weak  heterogeneity,  phase  distortion 
dominates  amplitude  distortion  [6]. 

We  propose  Adaptive  and  Robust  Techniques  (ART)  based 
on  RCB  for  TAT.  ART  can  be  used  to  mitigate  the  amplitude 
and  phase  distortion  problems  in  TAT  by  allowing  certain  un¬ 
certainties.  We  will  demonstrate  the  excellent  performance  of 
ART  by  using  both  data  simulated  on  a  2-D  breast  model  and 
data  experimentally  measured  from  mastectomy  specimens. 

II.  Problem  Formulation 

Consider  a  TAT  imaging  system  as  shown  in  Figure  1(a). 
A  stimulating  microwave  pulse  is  absorbed  by  the  biological 
tissue  under  testing,  which  causes  a  sudden  heat  change.  Due 
to  the  thermoacoustic  effect,  an  acoustic  pulse  is  generated 
which  can  be  recorded  by  an  ultrasonic  transducer  array  (a 
real  aperture  array  or  a  synthetic  aperture  array).  The  number 
of  transducers  in  the  array  (or  the  number  of  transducer  data 
acquisition  locations  in  the  synthetic  aperture  array)  is  M. 
Each  transducer  is  assumed  to  be  omnidirectional. 

The  data  model  for  the  sampled  and  digitized  acoustic  signal 
recorded  by  the  rath  transducer  is  given  by: 

xm(n)  =  sm(n)  +  em(n),  m  =  (1) 

where  n  is  the  discrete  time  index,  starting  from  t0  after 
excitation  pulse.  The  scalar  sm(n )  denotes  the  signal  com¬ 
ponent,  which  corresponds  to  the  acoustic  pulse  generated  at 
a  focal  point,  and  em(n)  is  the  residual  term,  which  includes 
unmodelled  noise  and  interference. 
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The  goal  of  ART  is  to  reconstruct  an  image  of  thermoa¬ 
coustic  response  intensity  /( r),  which  is  directly  related  to  the 
absorption  property  of  the  tissue,  from  the  recorded  data  set 
{xm(n)}.  Herein  the  vector  r  denotes  the  focal  point  location 
coordinate.  To  form  an  image,  we  scan  the  focal  point  at 
location  r  to  cover  the  entire  cross-section  of  the  tissue. 

The  discrete  arrival  time  of  the  pulse  (for  the  rath  trans¬ 
ducer)  can  be  determined  approximately  as: 


tm(r) 


±o_  llr  —  rm  || 
At  AtVQ 


(2) 


We  will  omit  the  dependence  of  the  arrival  time  tm( r)  on 
r  hereafter  for  notational  simplicity.  Here  At  is  the  sampling 
interval,  and  the  3-D  vector  rm  denotes  the  location  of  the  mth 
transducer.  The  sound  speed  v0  is  chosen  to  be  the  average 
sound  speed  of  the  biological  tissue  under  interrogation.  The 
notation  ||x||  denotes  the  Euclidean  norm  of  x,  and  \_y\  stands 
for  rounding  to  the  greatest  integer  less  than  y. 

The  signal  components  {sm(n)}^f=1  are  approximately 
scaled  and  shifted  versions  of  a  nominal  waveform  s(t)  at 
the  source: 


sm(n) 


exp(— a\\r  -  rm||) 
llr-rm|| 


s(n  -  tm), 


(3) 


where  a  is  the  attenuation  coefficient  in  Nepers/m.  We  prepro¬ 
cess  the  data  to  time  delay  all  the  signals  from  the  focal  point 
r  and  compensate  for  the  loss  in  amplitude  due  to  propagation 
decay.  Let  ym(n)  denote  the  signal  after  preprocessing  to 
backpropagate  the  detected  signal  to  the  source: 


Vm(n )  =  exp  (a||r  -  rm||)  •  ||r  -  rm|| xm(n  +  tm);  (4) 


Fig.  1.  (a):  A  schematic  of  a  2-D  synthetic  aperture  based  TAT  scanning 

system;  (b):  a  typical  acoustic  pulse  recorded  by  a  transducer  (for  data 
measured  from  breast  specimen  II). 


III.  Step  I  of  ART:  Waveform  Estimation 

The  first  step  of  ART  is  to  estimate  the  waveform  of  the 
acoustic  pulse  generated  by  the  focal  point  at  location  r,  based 
on  the  data  model  in  (5).  It  will  appear  that  we  have  neglected 
the  presence  of  phase  distortion  by  using  this  data  model  in 
the  first  step.  However,  by  allowing  a0  to  be  uncertain,  we  can 
tolerate  some  phase  distortions  as  well.  This  approximation 
causes  little  performance  degradation  to  our  robust  algorithm. 

Covariance  fitting  based  RCB  [4]  is  used  to  first  estimate 
the  steering  vector  ao,  and  use  the  estimated  ao  to  obtain 
an  optimal  beamformer  weight  vector  for  pulse  waveform 
estimation.  By  assuming  that  the  true  steering  vector  lies  in 
the  vicinity  of  the  nominal  steering  vector  a,  we  consider  the 
following  optimization  problem  [4]: 

max  <t2  subject  to  R  —  cr2aoaJ  >  0, 

cr2,a0 

||a0  -  a||2  <  e,  (6) 


then  the  received  vector  data  model  can  be  written  as: 

y(n)  =  aos(n)  +  e(n),  n  =  — N ,  •  •  • ,  N,  (5) 

where  ao  is  the  corresponding  steering  vector,  which 
is  approximately  equal  to  a  =  [1,  •  •  • ,  1]T,  y(n)  = 

[yi  (n),  •  •  • ,  Vm(^)]T?  e(n)  represents  the  noise  and  interfer¬ 
ence  term  after  preprocessing,  and  (-)T  denotes  the  transpose. 
The  time  interval  of  interests  for  the  signal  y (t)  is  from  —N 
to  N ,  which  is  large  enough  to  covers  the  expected  signal 
duration  in  the  region  of  interest. 

In  reality,  both  the  amplitude  and  the  phase  (or  pulse 
arrival  time)  of  the  acoustic  pulse  will  be  distorted.  Amplitude 
distortion  is  mainly  caused  by  multi-path  in  the  heterogeneous 
medium.  Phase  distortion  is  due  to  the  nonuniform  sound 
speed,  an  inaccurate  estimate  of  t0,  and  the  transducer  cal¬ 
ibration  error.  These  distortions  will  blur  the  image,  raise  the 
image  background  noise  level,  and  lower  the  values  of  the 
object  of  interest  [6]. 

In  ART,  we  mitigate  the  effects  of  these  distortions  by 
allowing  a0  to  belong  to  an  uncertainty  set  centered  at  a  and  by 
considering  the  signal  arriving  within  the  interval  from  —N  to 
N.  Our  ART  algorithm  consists  of  three  steps:  Step  I,  Robust 
Capon  Beamforming  (RCB)  [4]  for  robust  waveform  estima¬ 
tion;  Step  II,  peak-searching  for  phase  aberration  mitigation; 
Step  III,  intensity  calculation  for  forming  the  final  images. 


where  a2  is  the  power  of  the  signal  of  interest,  and 
1  N 

R  =  ^vTT  ^  y(n)yT(n)  E) 

is  the  sample  covariance  matrix.  The  second  constraint  in  (6) 
is  a  spherical  uncertainty  set;  an  elliptical  uncertainty  set  can 
be  used  instead  if  a  tighter  constraint  is  desirable  [4]. 

The  parameter  e  in  (6)  determines  the  size  of  the  uncertainty 
set  and  is  a  user  parameter.  To  avoid  the  trivial  solution  of 
a0  =  0,  we  require  that 

e<||a||2.  (8) 

To  attain  high  resolution  and  to  effectively  suppress  interfer¬ 
ence,  ^  should  be  made  as  small  as  possible.  On  the  other  hand, 
the  smaller  the  sample  size  N  or  the  larger  the  distortions,  the 
larger  should  5  be  chosen  [4].  Since  the  performance  of  RCB 
does  not  depend  very  critically  on  the  choice  of  5,  in  our 
examples  we  choose  certain  reasonable  initial  values  for  e, 
and  then  make  some  adjustments  empirically  based  on  image 
quality. 

By  using  the  Lagrange  multiplier  method,  the  solution  to 
(6)  is  given  by  [4]: 

a0  =  a  -  [)  +  fi r]  a,  (9) 
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where  I  is  the  identity  matrix,  fi  >  0  is  the  corresponding 
Lagrange  multiplier  that  can  be  solved  from  the  following 
equation: 

||(I  +  /xR)-1a||2  =e.  (10) 

Consider  the  eigendecomposition  on  the  sample  covariance 
matrix  R: 


r  =  uruT, 


(ii) 


where  the  columns  of  U  are  the  eigenvectors  of  R  and  the 
diagonal  matrix  T  consists  of  the  corresponding  eigenvalues 
7i  >  72  >  •  •  •  >  7m-  Let  b  m  UTa,  where  6m  denotes  its 
mth  element.  Then  (10)  can  be  rewritten  as 


M 


Am)  = 


\K 


^  (1  +  l-Hm ) ' 


= 


(12) 


Note  that  C(fi)  is  a  monotonically  decreasing  function  of  /i, 
with  £(0)  >  £  by  (8)  and  lim^oo  C(fi)  =  0  <  s,  which 
means  that  /jl  can  be  solved  efficiently,  say,  by  using  the 
Newton’s  method  (see  [4]  for  more  details).  After  obtaining 
the  value  of  /jl,  the  estimate  a0  of  the  actual  steering  vector  a0 
is  determined  by  (9). 

Observe  that  there  is  a  “scaling  ambiguity”  in  (6)  by 
treating  both  the  signal  power  a2  and  the  steering  vector  ao  as 
unknowns  (see  [4]).  To  eliminate  this  ambiguity,  we  scale  the 
solution  a0  to  make  its  norm  satisfy  the  following  condition: 


||ao||2  =  M. 


(13) 


To  obtain  an  estimate  for  the  signal  waveform  s(n),  we  ap¬ 
ply  a  weight  vector  to  the  preprocessed  signals  {y(n)}^=_N. 
The  weight  vector  is  determined  by  using  the  estimated 
steering  vector  a0  in  the  weight  vector  expression  of  the 
standard  Capon  beamformer  (see,  e.g.,  [4]): 

r  .  -i-l 


WRCB  = 


INI  _ 

M1/2  _7 


R  +  -I 


r+ n 


-1 


R  [r+  1I 

L  M  J 


-1 


(14) 


a0 


Note  that  (14)  has  a  diagonal  loading  form,  which  allows 
the  sample  covariance  matrix  to  be  rank-deficient.  The  beam- 
former  output  can  be  written  as: 


sRCB(n)  =  w^By(n),  n  =  -N,  ■  ■  ■ ,  N,  (15) 


which  is  the  waveform  estimate  for  the  acoustic  pulse  gener¬ 
ated  at  the  focal  point  at  location  r. 

RCB  can  provide  a  much  better  waveform  estimate  than  the 
conventional  DAS  but  at  a  higher  computational  cost.  For  a 
single  focal  point,  RCB  requires  0(M3)  flops,  which  mainly 
comes  from  the  eigen-decomposition  of  the  sample  covariance 
matrix  R  [4];  DAS  needs  only  O(M)  flops.  DAS  can  be 
used  as  a  fast  image  reconstruction  method  to  provide  initial 
imaging  results. 

The  weight  vector  used  by  DAS  for  waveform  estimation 
is 


wDAS  =  a, 


(16) 


and  the  estimated  waveform  is  given  by 

M 

sDAs(n)  =  w^sy (n)  =  ^  ym(n),  n  =  -N,  ■  ■  ■ ,  N.  (17) 

ra= 1 

IV.  Step  II  of  ART:  Peak  Searching 

Based  on  the  estimated  waveform  obtained  in  Step  I  for  the 
focal  point  at  location  r,  in  Step  II  of  ART,  we  will  search  for 
the  two  peaks  of  the  bipolar  acoustic  pulse  generated  by  the 
focal  point.  In  a  homogeneous  background  we  can  accurately 
calculate  the  arrival  time  of  the  acoustic  pulse  generated  by 
the  focal  point  at  location  r  by  using  (2);  this  is  not  true  in 
heterogeneous  biological  tissues  due  to  the  nonuniform  sound 
speed.  We  can  assume  that  the  original  peak  remains  a  peak 
in  the  waveform  estimated  from  Step  I  of  ART,  since  the 
amplitude  distortion  is  not  severe  in  TAT  [6]. 

The  bipolar  acoustic  pulse  has  one  peak  positive  and  another 
negative.  We  determine  the  positive  and  negative  peak  values 
as  follows: 

p+ =  max  <  max  s(n),0>,  (18) 

(ne[— A,A]  J 


P  =  min  <  min  s(n),0>,  (19) 

(nE[— A,A]  J 

where  the  searching  range  [—A,  A]  G  [—A,  TV]  is  around 
the  calculated  arrival  time  given  by  (2).  Here  A  is  a  user 
parameter.  Since  the  peak  searching  is  independent  of  the 
particular  waveform  estimation  methods,  we  use  s(n)  to 
denote  the  waveform  estimated  by  either  DAS  or  ART. 

The  search  range  is  determined  by  the  difference  between 
the  true  arrival  time  im  and  the  calculated  arrival  time  tm 
based  on  (2).  An  expression  the  arrival  time  difference  has 
been  given  in  [6]  in  the  following  form: 

5m(r')  =  tm  -  tm  oc  (20) 

where  r'  is  a  point  within  the  line  connecting  the  focal  point 
at  location  r  and  the  mth  transducer  at  location  rm,  and  u(r') 
is  the  local  sound  speed.  In  (20),  the  higher  order  terms  of 
[u(r')  —  vo]/ vq  have  been  ignored.  It  is  reasonable  to  assume 
that  u(r')  is  Gaussian  distributed  with  mean  v0  and  variance 
a2.  Consequently  the  arrival  time  difference  is  also  Gaussian 
distributed  with  zero-mean  and  variance  a2  oc  cf2/vq.  If  we 
choose  A  =  as,  and  the  duration  of  the  acoustic  pulse  is  r, 
we  can  find  the  two  peaks  of  the  pulse  within  the  interval 
(— as,  as  +  t)  on  the  recorded  signals  with  a  high  probability 
of  0.6826.  This  analysis  is  consistent  with  the  experimental 
measurements  in  [7].  From  our  examples,  we  found  that  a 
symmetric  range  [—A,  A]  around  the  estimated  arrival  time 
performs  similarly  to  the  asymmetric  range  [—A,  A  +  r],  and 
we  use  the  former  since  it  is  easy  to  handle  in  practice.  Also, 
we  can  use  similar  techniques  as  those  in  [7]  to  estimate  as 
to  find  a  good  searching  range  for  Step  II  of  ART,  and  to 
estimate  r  for  the  energy  type  methods  in  our  examples. 

There  is  a  tradeoff  in  choosing  the  searching  range.  The 
larger  the  searching  range,  the  higher  the  probability  we  can 
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find  the  peaks  of  the  acoustic  pulse  within  the  range.  However, 
if  the  range  is  chosen  too  large,  the  interferences  may  cause 
false  peaks,  and  as  a  consequence,  we  are  more  likely  to  find 
a  false  peak.  In  our  examples  we  choose  the  best  searching 
range  empirically  based  on  the  estimates  of  a#. 


V.  Step  III  of  ART:  Intensity  Calculation 


Based  on  the  estimation  of  the  waveform  generated  by  the 
focal  point  at  location  r,  the  response  intensity  for  need  to 
be  calculated.  There  are  two  major  types  of  response  intensity 
measurement  approaches:  amplitude  based  and  energy  based. 
They  extract  different  information  from  the  estimated  wave¬ 
form,  which  may  be  useful  to  physicians  in  different  ways. 

Conventional  DAS  uses  the  amplitude  based  measure  for 
TAT  imaging  [3],  with  the  corresponding  response  intensity 
given  by  5(0),  or  equivalently: 


M 

Ic  =  s( 0)  =  ^2  Vm{ 0). 

m=l 

The  energy  based  measure,  such  as  the  one  used 
mono-static  and  multi-static  microwave  imaging  for 
cancer  detection  [8],  calculate  the  intensity  by: 

r  t  ~  M  ^ 

^  =  j2(n)  =  J2  E  Vm(^>  > 

n— 0  n= 0  \_m=l 


(21) 

in  the 
breast 


(22) 


We  can  consider  using  the  peak  value  as  the  response 
intensity  measure  due  to  the  bipolar  nature  of  the  response 
at  the  focal  point: 


f  P+  if|P+|  >|P-|; 
|  P ~  otherwise, 


(23) 


Herein  we  keep  the  sign  of  the  maximum  amplitude  since  the 
sign  of  the  peak  may  also  contain  some  information  about  the 
focal  point. 

Peak-searching  maximizes  the  output  signal-to-noise  ratio. 
An  intuitive  explanation  is  that,  given  the  fact  that  the  acoustic 
pulse  is  bipolar  [3],  if  we  assume  that  the  residual  term  e(t)  is 
stationary,  or  its  power  is  uniform  over  time,  then  the  signal- 
to-noise  ratio  is  maximized  at  the  (positive  or  negative)  peak 
of  the  acoustic  pulse.  As  a  comparison,  the  conventional  DAS 
(21)  fixes  the  samples  to  be  summed  up  at  the  calculated  arrival 
time.  Due  to  phase  distortions,  the  waveform  at  the  calculated 
time  may  be  far  from  the  peak  value. 

We  can  also  employ  peak-to-peak  difference  as  the  response 
intensity  for  the  focal  point  at  location  r: 


Jpp  =  P+  -  P~  >  0.  (24) 


Peak-to-peak  difference  has  higher  imaging  contrast  than  peak 
value  measure:  the  peak-to-peak  difference  of  the  bipolar  pulse 
is  approximately  twice  the  absolute  peak  value,  which  means 
that  the  output  signal  power  of  the  former  is  four  times  of 
the  latter;  yet  the  noise  power  of  the  former  is  only  twice  of 
the  latter.  Therefore  the  output  Signal-to-Noise  Ratio  (SNR) 
is  doubled  by  using  the  peak-to-peak  difference  rather  than 
the  peak  value.  Both  peak-value  and  peak-to-peak  difference 
measures  belong  to  the  amplitude  based  measures. 


VI.  Numerical  and  Experimental  Examples 

A.  Numerical  Examples 

We  simulate  a  2-D  cylindrical  breast  model  using  the  Finite- 
Difference  Time-Domain  (FDTD)  method.  The  2-D  breast 
model  includes  2  mm  thick  skin,  chest  wall,  as  well  as 
randomly  distributed  fatty  breast  tissues  and  glandular  tissues. 
The  cross-section  of  the  breast  model  is  a  half  circle  with  a 
10  cm  diameter.  A  2  mm-diameter  tumor  is  located  at  2.2  cm 
below  the  skin  (at  x  =  7.0  cm,  y  =  6.0  cm). 

The  dielectric  properties  of  the  breast  tissues  are  assumed 
to  be  Gaussian  random  variables  with  variations  of  ±10% 
around  their  nominal  values.  The  nominal  values  are  chosen 
to  be  typical  of  those  reported  in  the  literature,  which  is 
given  in  Table  1  of  [8].  The  dispersive  properties  of  the 
tissues  are  also  considered  in  the  model.  Fat  breast  tissues  and 
glandular  tissues  are  randomly  distributed  in  the  breast  model. 
The  differential  equations  in  [9]  are  used  to  characterize  the 
thermoacoustic  effect  in  our  FDTD  simulations. 

We  model  the  sound  speed  within  the  breast  as  a  Gaussian 
random  variable  with  variation  ±5%  around  the  assumed 
average  sound  speed  of  1500  m/s.  The  probing  microwave 
pulse  used  here  is  a  modulated  rectangular  pulse  with  a 
modulating  frequency  of  800  MHz.  The  duration  of  the  pulse  is 
1  /as.  Before  applying  the  aforementioned  preprocessing  steps 
and  ART,  we  remove  the  strong  skin  response  using  techniques 
similar  to  those  in  [8].  In  the  following  all  the  images  are 
displayed  on  a  linear  scale,  and  we  will  name  the  imaging 
methods  by  their  waveform  estimation  method  followed  by 
the  intensity  calculation  approach,  such  as  “DAS-C.” 

Figure  2  shows  the  images  for  the  simulated  breast  model 
ART  and  DAS.  The  tumor  response  is  weak  for  such  a  small 
tumor.  In  these  images  we  use  e  =  0.1  M  and  the  searching 
range  [-14,  14].  Figure  2(a)  corresponds  to  DAS-C,  where 
the  tumor  is  buried  by  interference  and  noise.  In  Figure  2(b), 
for  RCB-E,  most  of  the  clutters  are  cleared  up  but  a  strong 
clutter  shows  up  near  the  chest  wall.  Images  produced  by 
ART-P  in  Figure  2(d)  and  by  ART-PP  in  Figure  2(f)  have 
lower  sidelobe  levels  and  higher  resolutions  than  their  DAS 
based  counterparts.  ART-PP  has  a  higher  contrast  than  the 
former,  due  to  the  peak-to-peak  difference  used  as  the  intensity 
measure. 

B.  Experimental  Results 

We  have  also  tested  ART  and  DAS  on  some  TAT  experi¬ 
mental  data  obtained  from  a  mastectomy  specimen  [2]  by  the 
Optical  Imaging  Faboratory  at  the  Texas  A&M  University.  In 
the  experiment,  the  breast  specimen  was  formed  to  a  cylin¬ 
drical  shape  inside  a  plastic  bowl.  The  bowl  was  immersed 
in  ultrasound  coupling  medium  in  a  container.  The  acoustic 
signals  were  recorded  at  240  equally  spaced  scanning  stops 
on  a  circular  track  of  radius  12.9  cm.  The  thickness  of  this 
specimen  was  about  4  cm  in  a  round  plastic  bowl  of  17  cm  in 
diameter.  This  lesion  was  diagnosed  as  an  invasive  metaplastic 
carcinoma  with  chondroid  and  squamous  metaplasia.  The  size 
of  the  tumor  was  measured  to  be  35  mm  in  diameter  by  TAT, 
and  36  mm  in  diameter  by  radiography  (see  [2]  for  details). 
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Fig.  2.  Reconstructed  images  based  on  the  2-D  simulated  breast  model  with 
one  2  mm-diameter  tumor. 


Figure  3  shows  the  reconstructed  images  for  the  breast 
specimen.  In  the  following  images,  the  searching  range  was  set 
to  [-3,  3],  and  e  =  0.5M  for  all  the  RCB  used  herein.  In  Figure 
3(a),  for  DAS-C,  the  dark  region  shows  an  blurred  object 
corresponding  to  the  breast  tumor.  In  3(b),  for  RCB-E,  the 
light  region  shows  a  vague  boundary  of  the  tumor.  In  Figures 
3(c),  for  DAS-P,  and  3(d),  for  ART-P,  a  dark  region  with  a 
clear  cut  has  a  good  correspondence  with  the  location  and 
shape  of  the  tumor  in  the  radiograph  3(g).  In  Figures  3(e),  for 
DAS-PP,  and  3(f),  for  ART-PP,  not  only  a  clear  image  of  the 
tumor  is  obtained,  but  also  the  detailed  boundary  is  revealed. 
For  comparison,  the  images  from  X-ray  mammography  and 
the  exact  inverse  solution  of  TAT  (see  [2]  for  more  details) 
are  shown  in  Figures  3(g)  and  3(h),  respectively. 

VII.  Conclusions 

Adaptive  and  Robust  Techniques  (ART)  have  been  proposed 
for  thermoacoustic  tomography.  ART  is  robust  to  the  am¬ 
plitude  and  phase  distortions  in  the  recorded  signals  caused 
by  the  acoustic  heterogeneity  of  biological  tissues.  Exam¬ 
ples  based  on  a  numerically  simulated  2-D  breast  model 
and  experimentally  measured  data  from  human  mastectomy 
specimens  demonstrate  the  excellent  performance  of  ART: 
high  resolution,  low-side  lobe  level,  much  better  interference 
suppression  capability.  ART  could  detect  tumors  as  small  as 
2-mm  in  diameter. 
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